data-analysis-and-classification-using-xgboost.ipynb

{"cells":[{"cell_type":"markdown","metadata":{"_cell_guid":"fa3d2a9c-442b-44c2-bf4d-45bb544ec369","_uuid":"bcba1675eebc008a35c21f6b64174161ff7d4c48"},"source":["# Sloan Digital Sky Survey Classification\n","## Classification of Galaxies, Stars and Quasars based on the RD14 from the SDSS"]},{"cell_type":"markdown","metadata":{"_cell_guid":"22caeac3-dc2a-4e42-a6c8-77d638f9e463","_uuid":"791df65ccbdbe25303b458851522f909ec27df61"},"source":["<h1 id=\"tocheading\">Table of Contents</h1>\n","<div id=\"toc\"></div>"]},{"cell_type":"code","execution_count":null,"metadata":{"_cell_guid":"09d820d4-f514-4162-8c14-6fb76f6ce985","_uuid":"fb94823b1affd20a31b1076e075fa86e9f106024","collapsed":true},"outputs":[],"source":["%%javascript\n","$.getScript('https://kmahelona.github.io/ipython_notebook_goodies/ipython_notebook_toc.js')"]},{"cell_type":"markdown","metadata":{"_cell_guid":"3179bd3c-9962-466d-89e4-ae7a9fce1893","_uuid":"f1b22e6594af55d511f16280928408646ac56462"},"source":["### About the notebook"]},{"cell_type":"markdown","metadata":{"_cell_guid":"d4bef474-9b3a-4bfd-80ee-b79f7fa144dc","_uuid":"73e8f01266e59393a11b6e00300eb916ee50662f"},"source":["In this notebook we will try to classify observations of space to be either stars, galaxies or quasars. We will try to have a complete cycle of the data science workflow including querying the database to get the dataset, data analysis and building machine learning models to predict for new data.\n","\n","We are using data from the Sloan Digital Sky Survey (**Release 14**). \n","\n","**I followed Niklas Donges' Titanic classification approach for this notebook as I found it really useful and comprehensive!**  "]},{"cell_type":"markdown","metadata":{"_cell_guid":"8bab71c2-754a-48ab-90cc-035977be89b8","_uuid":"56737fa09e0277f960d52b36c6589f08f799eeba"},"source":["### About the SDSS"]},{"cell_type":"markdown","metadata":{"_cell_guid":"3ae0f535-2034-41ce-8f96-ec53c9319738","_uuid":"710358dc8a224cd40398183cbb2081fcf8f96a50"},"source":["The Sloan Digital Sky Survey is a project which offers public data of space observations. Observations have been made since 1998 and have been made accessible to everyone who is interested. \n","\n","For this purpose a special 2.5 m diameter telescope was built at the Apache Point Observatory in New Mexico, USA. The telescope uses a camera of 30 CCD-Chips with 2048x2048 image points each. The chips are ordered in 5 rows with 6 chips in each row. Each row observes the space through different optical filters (u, g, r, i, z) at wavelengths of approximately 354, 476, 628, 769, 925 nm.\n","\n","The telescope covers around one quarter of the earth's sky - therefore focuses on the northern part of the sky.\n","\n","**For more information about this awesome project - please visit their website:**\n","\n","http://www.sdss.org/\n","\n","![alt text](http://www.fingerprintdigitalmedia.com/wp-content/uploads/2014/08/sdss1.jpg)"]},{"cell_type":"markdown","metadata":{"_cell_guid":"b9d4b609-9f0f-4296-8462-98975b2ece09","_uuid":"16d228ddb3b0d71d6e13093552c04a21146b75e5"},"source":["### Importing Libraries"]},{"cell_type":"code","execution_count":null,"metadata":{"_cell_guid":"13914766-c2fb-4801-8846-6c78e6d1cb03","_uuid":"5bb212bdb34abc34f8bed1a0bc2d1a6287166221","collapsed":true},"outputs":[],"source":["import numpy as np\n","import pandas as pd\n","import matplotlib.pyplot as plt\n","from matplotlib import style\n","import seaborn as sns\n","sns.set_style('whitegrid')\n","import tensorflow as tf\n","from sklearn.model_selection import train_test_split, cross_val_predict\n","from sklearn.svm import SVC\n","from sklearn.ensemble import RandomForestClassifier\n","from sklearn.linear_model import SGDClassifier\n","from sklearn.neighbors import KNeighborsClassifier\n","from sklearn.naive_bayes import GaussianNB\n","from sklearn.decomposition import PCA\n","from xgboost import XGBClassifier\n","from sklearn.preprocessing import LabelEncoder, MinMaxScaler\n","from sklearn.metrics import confusion_matrix, precision_score, recall_score, f1_score\n","import time\n","import warnings\n","warnings.filterwarnings(\"ignore\", category=DeprecationWarning)\n","%matplotlib inline\n","\n","SMALL_SIZE = 10\n","MEDIUM_SIZE = 12\n","\n","plt.rc('font', size=SMALL_SIZE)\n","plt.rc('axes', titlesize=MEDIUM_SIZE)\n","plt.rc('axes', labelsize=MEDIUM_SIZE)\n","plt.rcParams['figure.dpi']=150"]},{"cell_type":"markdown","metadata":{"_cell_guid":"87680814-020a-4a48-966a-188ffb393eb8","_uuid":"70ed9474476ba00cbe256cfd857e662518ba5d14"},"source":["## Data Acquisition"]},{"cell_type":"markdown","metadata":{"_cell_guid":"c9336c15-46c3-4a2b-8ce3-002f490aad7c","_uuid":"32db81f80cb49de30a8a6c3ad04972d9de971c0a"},"source":["Public data from the SDSS can be accessed through multiple ways - I used the **CasJobs** website which offers a **SQL-based interface** which lets you query their database which contains the released data.\n","\n","For more information about how to get data from the SDSS see their Data Access Guide:\n","\n","http://www.sdss.org/dr14/data_access/\n","\n","I used the sample query given by the **CasJobs** to receive the data. Find the exact query below:"]},{"cell_type":"markdown","metadata":{"_cell_guid":"30913682-a98d-4910-a9a1-e2b2cbe01c9e","_uuid":"be5d31d91aa2035ef2a3198e1a2fc1165fe66da2"},"source":["### Query"]},{"cell_type":"markdown","metadata":{"_cell_guid":"34d34aba-b697-4114-bc08-68d9374325df","_uuid":"e50d82fc8cca0b72bca172cf4848b1c2cfc0d8a3"},"source":["**SELECT TOP 10000** <br/>\n","p.objid,p.ra,p.dec,p.u,p.g,p.r,p.i,p.z, p.run, p.rerun, p.camcol, p.field,  <br/>\n","s.specobjid, s.class, s.z as redshift, s.plate, s.mjd, s.fiberid  <br/>\n","**FROM** PhotoObj **AS** p <br/>\n","   **JOIN** SpecObj **AS** s **ON** s.bestobjid = p.objid <br/>\n","**WHERE** <br/>\n","   p.u **BETWEEN** 0 **AND** 19.6 <br/>\n","   **AND** g **BETWEEN** 0 **AND** 20 <br/>\n","\n"]},{"cell_type":"markdown","metadata":{"_cell_guid":"a789e0cf-25fa-440c-bc40-cc1b63a571da","_uuid":"ff831f4235020971246b431da8e86fb2fa2a8986"},"source":["The above query joins two tables (actually views): The image table (PhotoObj) which contains all image objects and the spectral table (SpecObj) which contains corresponding spectral data. "]},{"cell_type":"markdown","metadata":{"_cell_guid":"47aefcd7-8ce4-41b8-848c-610ec48aa89c","_uuid":"2ad8534b4a25b74d0c0e162949202b1e9de206a9"},"source":["### Feature Description"]},{"cell_type":"markdown","metadata":{"_cell_guid":"2b272b59-8c37-4e8e-a78a-ac5f50919d2e","_uuid":"82c44461c45ae463a3d429f564f8bca17f8b0f4c"},"source":["#### View \"PhotoObj\"\n","* objid = Object Identifier\n","* ra = J2000 Right Ascension (r-band)\n","* dec = J2000 Declination (r-band)\n","\n","Right ascension (abbreviated RA) is the angular distance measured eastward along the celestial equator from the Sun at the March equinox to the hour circle of the point above the earth in question. When paired with declination (abbreviated dec), these astronomical coordinates specify the direction of a point on the celestial sphere (traditionally called in English the skies or the sky) in the equatorial coordinate system.\n","\n","Source: https://en.wikipedia.org/wiki/Right_ascension\n","\n","* u = better of DeV/Exp magnitude fit\n","* g = better of DeV/Exp magnitude fit\n","* r = better of DeV/Exp magnitude fit\n","* i = better of DeV/Exp magnitude fit\n","* z = better of DeV/Exp magnitude fit\n","\n","The Thuan-Gunn astronomic magnitude system. u, g, r, i, z represent the response of the 5 bands of the telescope.\n","\n","Further education: https://www.astro.umd.edu/~ssm/ASTR620/mags.html\n","\n","* run = Run Number\n","* rereun = Rerun Number\n","* camcol = Camera column\n","* field = Field number\n","\n","Run, rerun, camcol and field are features which describe a field within an image taken by the SDSS. A field is basically a part of the entire image corresponding to 2048 by 1489 pixels. A field can be identified by:\n","- **run** number, which identifies the specific scan,\n","- the camera column, or \"**camcol**,\" a number from 1 to 6, identifying the scanline within the run, and\n","- the **field** number. The field number typically starts at 11 (after an initial rampup time), and can be as large as 800 for particularly long runs.\n","- An additional number, **rerun**, specifies how the image was processed. \n","\n","#### View \"SpecObj\"\n","\n","* specobjid = Object Identifier\n","* class = object class (galaxy, star or quasar object)\n","\n","The class identifies an object to be either a galaxy, star or quasar. This will be the response variable which we will be trying to predict.\n","\n","* redshift = Final Redshift\n","* plate = plate number\n","* mjd = MJD of observation\n","* fiberid = fiber ID\n","\n","In physics, **redshift** happens when light or other electromagnetic radiation from an object is increased in wavelength, or shifted to the red end of the spectrum. \n","\n","Each spectroscopic exposure employs a large, thin, circular metal **plate** that positions optical fibers via holes drilled at the locations of the images in the telescope focal plane. These fibers then feed into the spectrographs. Each plate has a unique serial number, which is called plate in views such as SpecObj in the CAS.\n","\n","**Modified Julian Date**, used to indicate the date that a given piece of SDSS data (image or spectrum) was taken.\n","\n","The SDSS spectrograph uses optical fibers to direct the light at the focal plane from individual objects to the slithead. Each object is assigned a corresponding **fiberID**. \n","\n","**Further information on SDSS images and their attributes:** \n","\n","http://www.sdss3.org/dr9/imaging/imaging_basics.php\n","\n","http://www.sdss3.org/dr8/glossary.php"]},{"cell_type":"code","execution_count":null,"metadata":{"_cell_guid":"ffb06ef6-73f7-4f42-ab42-5d5b5f773ba7","_uuid":"04e88f8c9c12167a1c23e47b3e2046246510e983","collapsed":true},"outputs":[],"source":["#sdss_df = pd.read_csv('Skyserver_SQL2_27_2018 6_51_39 PM.csv', skiprows=1)\n","sdss_df = pd.read_csv('../input/Skyserver_SQL2_27_2018 6_51_39 PM.csv', skiprows=1)"]},{"cell_type":"markdown","metadata":{"_cell_guid":"ca3f2be5-d15c-47ef-bc29-9a254ce2199d","_uuid":"776947da30d7e7e2b88bd87dca65eaa8b04b158c"},"source":["## Data Exploration"]},{"cell_type":"markdown","metadata":{"_cell_guid":"c5f29b18-e5f3-4b99-8abe-ce6f492f0584","_uuid":"8fd7a1ae6fb17276694385e527fbfd98d548be4a"},"source":["### Basic stats about our dataset"]},{"cell_type":"markdown","metadata":{"_cell_guid":"6246cf67-967f-4e01-99d0-191b0ea3a8d9","_uuid":"5f49f9775c7d5cb9e46e725dbceb1269d7c23416"},"source":["Let's take a first look at our dataset to see what we're working with!"]},{"cell_type":"code","execution_count":null,"metadata":{"_cell_guid":"26b460b4-d472-440b-b690-d8c12f267944","_uuid":"149d00bab55383b7014daf916df81410717b643a","collapsed":true,"scrolled":true},"outputs":[],"source":["sdss_df.head()"]},{"cell_type":"markdown","metadata":{"_cell_guid":"aec349a7-7dcd-4e93-b9d1-5a8f032ca6c7","_uuid":"43d8f5997c96f67ff4445ba00cd60c1d3f36d744"},"source":["We can tell that we have all the features as described in the above query. \n","\n","We notice that there are no categorical features at all - besides the class column. As some machine learning models can't handle categorical feature columns at all, we will encode this column to be a numerical column later on."]},{"cell_type":"markdown","metadata":{"_cell_guid":"811c152a-331b-4b69-a7d3-ab52934ef15c","_uuid":"7b054433b30b6ea6b1f9024ef323270d86874e76"},"source":["Let's find out about the types of columns we have:"]},{"cell_type":"code","execution_count":null,"metadata":{"_cell_guid":"bd29718c-d53e-4b52-a9ee-b8a9a252d864","_uuid":"926e1713524a5fcac0353e1d7c03d1f74042f348","collapsed":true,"scrolled":false},"outputs":[],"source":["sdss_df.info()"]},{"cell_type":"markdown","metadata":{"_cell_guid":"75a847f1-0dfc-4228-9cbc-49d6034463e5","_uuid":"9b01bc847e158cfa00d411ea687cb573a0037aef"},"source":["The dataset has 10000 examples, 17 feature columns and 1 target column. 8 of the 17 features are 64 bit integers, 1 feature is an unsigned 64 bit integer, 8 are 64 bit floats and the target column is of the type object. "]},{"cell_type":"code","execution_count":null,"metadata":{"_cell_guid":"6723c745-8446-46f0-a866-8c22668607d3","_uuid":"77495f8526975b41e2ba43063b82d807e8ba1109","collapsed":true,"scrolled":true},"outputs":[],"source":["sdss_df.describe()"]},{"cell_type":"markdown","metadata":{"_cell_guid":"8c6d5703-4e91-4c13-9a88-235c9848f6e8","_uuid":"6da0473e16af044373072b21e28d3425f6541cca"},"source":["From the above table we can tell that are no missing values at all. This means: **no imputing!**\n","\n","We also notice that most of the features stay within a reasonable scale when comparing values within **only one** column. We can recognize this from the min, max and quartil rows."]},{"cell_type":"code","execution_count":null,"metadata":{"_cell_guid":"980baf01-358b-4f40-b1cd-886a8463befe","_uuid":"75724e4e2a11a567b6155503ce0671ae3f5e2e7d","collapsed":true,"scrolled":true},"outputs":[],"source":["sdss_df['class'].value_counts()"]},{"cell_type":"markdown","metadata":{"_cell_guid":"47229348-1a06-4346-949b-1b2b135a20e9","_uuid":"d03d6a82a73b4889fc9bfee64e0773acd9b09302"},"source":["The most objects (50%) are galaxies, a little less (40%) are stars and only around (10%) of the rows are classified as QSOs."]},{"cell_type":"markdown","metadata":{"_cell_guid":"f890c211-edfb-45d5-be7d-545eb358220b","_uuid":"529e41f022c84c223fc508224b63b527677a61fc"},"source":["### First Data Filtering"]},{"cell_type":"code","execution_count":null,"metadata":{"_cell_guid":"31099fbb-fc0a-4956-bb26-dde90eecbb82","_uuid":"cb352ba3c067c076162c8ce399d8bcd4fe24f9c6","collapsed":true},"outputs":[],"source":["sdss_df.columns.values"]},{"cell_type":"markdown","metadata":{"_cell_guid":"25d1cf21-621e-462b-8d8f-faae4daca4d8","_uuid":"4eed688c6d406af4e335c9237431813b461fb636"},"source":["There is no need to know everything about stars, galaxy or quasars - yet we can already tell which features are **unlikely** to be related to the target variable 'class'.\n","\n","**objid** and **specobjid** are just identifiers for accessing the rows back when they were stored in the original databank. Therefore we will not need them for classification as they are not related to the outcome.\n","\n","Even more: The features 'run', 'rerun', 'camcol' and 'field' are values which describe parts of the camera at the moment when making the observation, e.g. 'run' represents the corresponding scan which captured the oject.\n","\n","Source: http://www.sdss3.org/dr9/imaging/imaging_basics.php\n","\n","We will drop these columns as any correlation to the outcome would be coincidentally."]},{"cell_type":"code","execution_count":null,"metadata":{"_cell_guid":"d1246406-e5f7-4cb7-8365-50c8e44e19d9","_uuid":"4d9519ab0d8b9f843219ce89531c3fc13dafc5be","collapsed":true,"scrolled":false},"outputs":[],"source":["sdss_df.drop(['objid', 'run', 'rerun', 'camcol', 'field', 'specobjid'], axis=1, inplace=True)\n","sdss_df.head(1)"]},{"cell_type":"markdown","metadata":{"_cell_guid":"6b86afe2-f926-416a-afd9-8b36edf9bbdd","_uuid":"bc4bc3db7af9a1e46c21c5fda26e023137e0e239"},"source":["### Univariate Analysis"]},{"cell_type":"markdown","metadata":{"_cell_guid":"09cc084c-d589-4393-aa72-3c60c60e2225","_uuid":"737c3fc575497bd4c675c4cf5949edfd56f87fd7"},"source":["#### Redshift"]},{"cell_type":"markdown","metadata":{"_cell_guid":"2fbe3492-24e7-4368-9382-9a0567ff9b09","_uuid":"3e55f4fbff1621e5d562e9353c646b447a1e003b"},"source":["To start the univariate analysis we will plot histograms for the 'redshift' feature column for each class.\n","\n","This will tell us how the redshift values are distributed over their range."]},{"cell_type":"code","execution_count":null,"metadata":{"_cell_guid":"f7e55b21-caa1-4b91-8f5a-c465d56cc146","_uuid":"5a6a061a52d8e4bb3286be05e8b31c0c656f3dfd","collapsed":true,"scrolled":true},"outputs":[],"source":["fig, axes = plt.subplots(nrows=1, ncols=3,figsize=(16, 4))\n","ax = sns.distplot(sdss_df[sdss_df['class']=='STAR'].redshift, bins = 30, ax = axes[0], kde = False)\n","ax.set_title('Star')\n","ax = sns.distplot(sdss_df[sdss_df['class']=='GALAXY'].redshift, bins = 30, ax = axes[1], kde = False)\n","ax.set_title('Galaxy')\n","ax = sns.distplot(sdss_df[sdss_df['class']=='QSO'].redshift, bins = 30, ax = axes[2], kde = False)\n","ax = ax.set_title('QSO')"]},{"cell_type":"markdown","metadata":{"_cell_guid":"a44d0c69-fb05-4b09-99c3-d92dc1ec0b9a","_uuid":"f389689f8d958bdd4e281089036af6d26c1c4d10"},"source":["This is an interesting result.\n","\n","We can cleary tell that the redshift values for the classes quite differ. \n","\n","* **Star:** The histogram looks like a truncated zero-centered normal distribution.\n","\n","* **Galaxy:** The redshift values may come from a slightly right-shifted normal distribution which is centered around 0.075.\n","\n","* **QSO:** The redshift values for QSOs are a lot more uniformly distributed than for Stars or Galaxies. They are roughly evenly distributed from 0 to 3, than the occurences decrease drastically. For 4 oder ~5.5 there are some outliers.\n","\n","**The redshift can be an estimate(!) for the distance from the earth to a object in space.**\n","\n","Hence the distplot tells us that most of the stars observed are somewhat closer to the earth than galaxies or quasars. Galaxies tend to be a little further away and quasars are distant from very close to very far.  \n","\n","Possible rookie explanation: Since galaxies and quasars radiate stronger due to their size and physical structure, they can be observed from further away than \"small\" stars.\n","\n","As we can distinct the classes from each other just based on this column - 'redshift' is very likely to be helping a lot classifying new objects."]},{"cell_type":"markdown","metadata":{"_cell_guid":"b6a6c486-ca64-47cb-ad60-298338949c32","_uuid":"76749b9ea7daeba52d051df437f67c723b5231e7"},"source":["#### dec"]},{"cell_type":"markdown","metadata":{"_cell_guid":"4fe6d3c7-336f-4ab7-95fb-f0c283f2e6a6","_uuid":"7bef264989ac504923ba4606da2c4d4aa157f59d"},"source":["Let's lvplot the values of dec (Recall: position on celestial equator)!"]},{"cell_type":"code","execution_count":null,"metadata":{"_cell_guid":"267fdf63-0a7c-4e8b-895a-3358d17202b5","_uuid":"e5eadada1e83442db6e878c79d12d798ed163d65","collapsed":true},"outputs":[],"source":["fig, axes = plt.subplots(nrows=1, ncols=1,figsize=(16, 4))\n","ax = sns.lvplot(x=sdss_df['class'], y=sdss_df['dec'], palette='coolwarm')\n","ax.set_title('dec')"]},{"cell_type":"markdown","metadata":{"_cell_guid":"9b580f6c-c6eb-45af-b741-56952e230e89","_uuid":"582f48aae4db70793a6f58247ee060b4edc7d6b0"},"source":["**First of all: what does this plot tell us?**\n","\n","The Letter value (LV) Plot show us an estimate of the distribution of the data. It shows boxes which relate to the amount of values within the range of values inside the box.\n","\n","In this case we can observe a clear distinction between Stars and the other two classes. The difference between Galaxies and Quasars is smaller.\n","\n","* **Star:** The largest part of the data points lay within a 0 to 10 range. Another large part consists of values between about 10 to 55. Only small amounts of the data are lower or higher than these ranges.\n","\n","* **Galaxy:** The largest part of values lays between 0 and 45. There is a smaller amount of values in the range of 45 to 60. The rest of the data has smaller or higher values.\n","\n","* **QSO:** This plot looks quite similiar to the GALAXY plot. Only the amount of data points in the range of 0 to 60 is even bigger.\n","\n","Side Note: The fact that the distribution of dec values of galaxies und quasar objects is almost the same might indicate that one can find both galaxies and quasars at smiliar positions in the night sky."]},{"cell_type":"markdown","metadata":{"_cell_guid":"83b72276-7e91-40bd-b450-cb5f05b2a5e9","_uuid":"fd26a7feb55ad882952611adced5aa71de01ad05"},"source":["### Multivariate Analysis"]},{"cell_type":"markdown","metadata":{"_cell_guid":"f3174f32-1e96-47d5-ab74-2a14f1e7b1f7","_uuid":"32ac7f6965966d8f25013c5a2e1d2817d9e78cf3"},"source":["#### u,g,r,i,z filters"]},{"cell_type":"markdown","metadata":{"_cell_guid":"27c33528-10f3-400f-a59a-b0c28416e00a","_uuid":"4761d9058cb828d6e016fe43693019981ce22596"},"source":["Recall: u, g, r, i, z represent the different wavelengths which are used to capture the observations.\n","\n","Let's find out how much they are correlated."]},{"cell_type":"code","execution_count":null,"metadata":{"_cell_guid":"4a856729-40ba-4af5-800a-855e7d866ae6","_uuid":"3b0a813154788f955256e31ea6984ca44074ac88","collapsed":true,"scrolled":true},"outputs":[],"source":["fig, axes = plt.subplots(nrows=1, ncols=3,figsize=(16, 4))\n","fig.set_dpi(100)\n","ax = sns.heatmap(sdss_df[sdss_df['class']=='STAR'][['u', 'g', 'r', 'i', 'z']].corr(), ax = axes[0], cmap='coolwarm')\n","ax.set_title('Star')\n","ax = sns.heatmap(sdss_df[sdss_df['class']=='GALAXY'][['u', 'g', 'r', 'i', 'z']].corr(), ax = axes[1], cmap='coolwarm')\n","ax.set_title('Galaxy')\n","ax = sns.heatmap(sdss_df[sdss_df['class']=='QSO'][['u', 'g', 'r', 'i', 'z']].corr(), ax = axes[2], cmap='coolwarm')\n","ax = ax.set_title('QSO')"]},{"cell_type":"markdown","metadata":{"_cell_guid":"b9f174fd-eca7-4a1c-a412-e58a2da4bf3d","_uuid":"0317928d84b4f79e74e611ceb07a7a4295c05fe5"},"source":["Right of the top we observe that the correlation matrices look very similiar for every class.\n","\n","We can tell that there are high correlations between the different bands. This feels not really suprising - intuitively one would think that if one of the bands captures some object, the other bands should capture something aswell.\n","\n","Therefore it is interesting to see that band 'u' is less correlated to the other bands. \n","\n","Remember: u, g, r, i, z capture light at wavelengths of 354, 476, 628, 769 and 925 nm.\n","\n","This might indicates that galaxies, stars and quasar objects shine brighter at wavelengths from 476 - 925 nm. Don't quote me on that though.\n","\n","**But:** as we can see - the correlation is roughly the same for every class...the different bands behave the same for the different classes!"]},{"cell_type":"markdown","metadata":{"_cell_guid":"0cbe1c99-0148-4406-a11d-17d1efa44025","_uuid":"7b2155d1a1261d2065a92008cb2a33d425f7b3dd"},"source":["#### Right ascension (ra) and declination (dec) "]},{"cell_type":"markdown","metadata":{"_cell_guid":"e4dac296-e9a1-4a63-83fc-46c3ca2265a2","_uuid":"4e248ae503d23d0dd3f2473108007ad7b899b47a"},"source":["We will now plot the right ascension versus the declination depending on the class "]},{"cell_type":"code","execution_count":null,"metadata":{"_cell_guid":"66068c6e-e97b-4e00-8baf-85d67ee6b155","_uuid":"6e0aa853466d5c45afec6ac1efd4f9db2a013642","collapsed":true,"scrolled":false},"outputs":[],"source":["sns.lmplot(x='ra', y='dec', data=sdss_df, hue='class', fit_reg=False, palette='coolwarm', size=6, aspect=2)\n","plt.title('Equatorial coordinates')"]},{"cell_type":"markdown","metadata":{"_cell_guid":"c61cc9f2-06eb-460e-a44e-84aafb16a3fa","_uuid":"3abe9b02be5e96d72dd5f08a5708fdaa55fd4775"},"source":["As we can clearly observe the equatorial coordinates do not differ significantly between the 3 classes. There are some outliers for stars and galaxies but for the bigger part the coordinates are within the same range.\n","\n","Why is that?\n","\n","All SDSS images cover the same area of the sky. The plot above tells us that stars, galaxies and quasars are observed equally at all coordinates within this area. So whereever the SDSS \"looks\" - the chance of observing a star or galaxy or quasar is always the same.  \n","\n","**This contradicts our interpretation of the letter value plot of dec from the univariate analysis.**"]},{"cell_type":"markdown","metadata":{"_cell_guid":"aa0a41b4-c96c-4193-8b86-fee485af4361","_uuid":"4e3e3eccee7ba9f0f2bcb5b562bda8b4d3c5a424"},"source":["## Feature Engineering"]},{"cell_type":"markdown","metadata":{"_cell_guid":"cf9c9fc8-d5cc-48c5-b79f-2c0212e6bb44","_uuid":"355a3ba3922a9553c8e79fc9309f7dc03b3a51d8"},"source":["### u, g, r, i, z"]},{"cell_type":"markdown","metadata":{"_cell_guid":"10fc6f0d-fad0-44ca-8996-bdf513e48358","_uuid":"2fe1eb4a169a3400d65c3f5d3134ae8fb250240e"},"source":["We will now reduce the amount of dimensions by replacing the different bands 'u', 'g', 'r', 'i' and 'z' by a linear combination with only 3 dimensions using **Principal Component Analysis**.\n","\n","**Principal Component Analysis:**\n","\n","n observations with p features can be interpreted as n points in a p-dimensional space. PCA aims to project this space into a q-dimensional subspace (with q<p) with as little information loss as possible. \n","\n","It does so by finding the q directions in which the n points vary the most (--> the principal components). It then projects the original data points into the q-dimensional subspace. PCA returns a n x q dimensional matrix. \n","\n","Using PCA on our data will decrease the amount of operations during training and testing."]},{"cell_type":"code","execution_count":null,"metadata":{"_cell_guid":"fbb45762-272f-40b8-bde6-9d3dd9c1cd55","_uuid":"8a97dca248a7b0473c784af669ea00b59017fa8a","collapsed":true},"outputs":[],"source":["sdss_df_fe = sdss_df\n","\n","# encode class labels to integers\n","le = LabelEncoder()\n","y_encoded = le.fit_transform(sdss_df_fe['class'])\n","sdss_df_fe['class'] = y_encoded\n","\n","# Principal Component Analysis\n","pca = PCA(n_components=3)\n","ugriz = pca.fit_transform(sdss_df_fe[['u', 'g', 'r', 'i', 'z']])\n","\n","# update dataframe \n","sdss_df_fe = pd.concat((sdss_df_fe, pd.DataFrame(ugriz)), axis=1)\n","sdss_df_fe.rename({0: 'PCA_1', 1: 'PCA_2', 2: 'PCA_3'}, axis=1, inplace = True)\n","sdss_df_fe.drop(['u', 'g', 'r', 'i', 'z'], axis=1, inplace=True)\n","sdss_df_fe.head()"]},{"cell_type":"markdown","metadata":{"_cell_guid":"4e797850-cbb6-44f8-84d8-2649d9dcb461","_uuid":"73423e72a3e388e8d6b8f16672a7b77215359cfe"},"source":["## Machine Learning Models - Training"]},{"cell_type":"markdown","metadata":{"_cell_guid":"64793a20-a257-4717-a53b-ab02c6cdffb1","_uuid":"e970946f57bff18edaed172d867bec65846f7413"},"source":["#### Feature Scaling"]},{"cell_type":"markdown","metadata":{"_cell_guid":"0a8fc626-e05b-483e-bed1-17c5be108903","_uuid":"282a92a1d893a96b40b70b1c6d1fae063ba08ba3"},"source":["We will now train different models on this dataset. \n","\n","Scaling all values to be within the (0, 1) interval will reduce the distortion due to exceptionally high values and make some algorithms converge faster."]},{"cell_type":"code","execution_count":null,"metadata":{"_cell_guid":"23fc8398-331d-4164-8779-0516264ece29","_uuid":"c147b9e899cc2dd508d736274c856f88fb49321d","collapsed":true},"outputs":[],"source":["scaler = MinMaxScaler()\n","sdss = scaler.fit_transform(sdss_df_fe)"]},{"cell_type":"markdown","metadata":{"_cell_guid":"13e1fd6f-820c-4098-a486-0079b300e6c9","_uuid":"e76e964a4b93c27e1ab01c24a06be09d8b304970"},"source":["We will  split the data into a training and a test part. The models will be trained on the training data set and tested on the test data set"]},{"cell_type":"code","execution_count":null,"metadata":{"_cell_guid":"a94f2305-4f2a-4d11-b368-528820d28b0b","_uuid":"da9a96996b695f8ca2c268f43106e179ba6be1f5","collapsed":true},"outputs":[],"source":["X_train, X_test, y_train, y_test = train_test_split(sdss_df_fe.drop('class', axis=1), \n","                                                    sdss_df_fe['class'], test_size=0.33)"]},{"cell_type":"markdown","metadata":{"_cell_guid":"d74dbe6b-ebf8-4b61-b568-c2dabc031e5c","_uuid":"3d19792d9ae5a949d2c524fb0c701a4084051816"},"source":["#### K Nearest Neighbors"]},{"cell_type":"code","execution_count":null,"metadata":{"_cell_guid":"a16b9525-b3c9-46cf-8714-12707caf2502","_uuid":"817a9fd86ca3d5195046fed7a72750538dee3051","collapsed":true},"outputs":[],"source":["knn = KNeighborsClassifier()\n","training_start = time.perf_counter()\n","knn.fit(X_train, y_train)\n","training_end = time.perf_counter()\n","prediction_start = time.perf_counter()\n","preds = knn.predict(X_test)\n","prediction_end = time.perf_counter()\n","acc_knn = (preds == y_test).sum().astype(float) / len(preds)*100\n","knn_train_time = training_end-training_start\n","knn_prediction_time = prediction_end-prediction_start\n","print(\"Scikit-Learn's K Nearest Neighbors Classifier's prediction accuracy is: %3.2f\" % (acc_knn))\n","print(\"Time consumed for training: %4.3f seconds\" % (knn_train_time))\n","print(\"Time consumed for prediction: %6.5f seconds\" % (knn_prediction_time))"]},{"cell_type":"markdown","metadata":{"_cell_guid":"cd50994d-fa25-4866-8a61-9d316a763a4a","_uuid":"e901dbdaed422f26f42d9396a3b7a668c876f401"},"source":["#### Naive Bayes"]},{"cell_type":"code","execution_count":null,"metadata":{"_cell_guid":"2db753f8-9832-4b47-9ba8-35db1448602f","_uuid":"39aa691fa46edc353cf096543a19fa9d9bc1498f","collapsed":true},"outputs":[],"source":["gnb = GaussianNB()\n","training_start = time.perf_counter()\n","gnb.fit(X_train, y_train)\n","training_end = time.perf_counter()\n","prediction_start = time.perf_counter()\n","preds = gnb.predict(X_test)\n","prediction_end = time.perf_counter()\n","acc_gnb = (preds == y_test).sum().astype(float) / len(preds)*100\n","gnb_train_time = training_end-training_start\n","gnb_prediction_time = prediction_end-prediction_start\n","print(\"Scikit-Learn's Gaussian Naive Bayes Classifier's prediction accuracy is: %3.2f\" % (acc_gnb))\n","print(\"Time consumed for training: %4.3f seconds\" % (gnb_train_time))\n","print(\"Time consumed for prediction: %6.5f seconds\" % (gnb_prediction_time))"]},{"cell_type":"markdown","metadata":{"_cell_guid":"026ded4a-6c4d-45ed-819e-68fe0f3e7238","_uuid":"e13c5a84383fc5e96ae786757f68a6d5d4fefb83"},"source":["#### XGBoost"]},{"cell_type":"code","execution_count":null,"metadata":{"_cell_guid":"ebc9caed-0883-4352-9fd8-dc40154416f3","_uuid":"e6f31bfdb40799c40b90ee4b603e6ffae9197f4d","collapsed":true},"outputs":[],"source":["xgb = XGBClassifier(n_estimators=100)\n","training_start = time.perf_counter()\n","xgb.fit(X_train, y_train)\n","training_end = time.perf_counter()\n","prediction_start = time.perf_counter()\n","preds = xgb.predict(X_test)\n","prediction_end = time.perf_counter()\n","acc_xgb = (preds == y_test).sum().astype(float) / len(preds)*100\n","xgb_train_time = training_end-training_start\n","xgb_prediction_time = prediction_end-prediction_start\n","print(\"XGBoost's prediction accuracy is: %3.2f\" % (acc_xgb))\n","print(\"Time consumed for training: %4.3f\" % (xgb_train_time))\n","print(\"Time consumed for prediction: %6.5f seconds\" % (xgb_prediction_time))"]},{"cell_type":"markdown","metadata":{"_cell_guid":"75b6b173-0875-4716-b458-7b741295fcf2","_uuid":"d28d0e8f0d89945c4cb989b19d1829a50987ccdc"},"source":["#### Scitkit-Learn's Random Forest Classifier"]},{"cell_type":"code","execution_count":null,"metadata":{"_cell_guid":"b5d4a318-08c2-4172-8755-a82aa4e0cafc","_uuid":"9dc11189c6289f7b3f0bb982bdf0c6b82da3e8d3","collapsed":true},"outputs":[],"source":["rfc = RandomForestClassifier(n_estimators=10)\n","training_start = time.perf_counter()\n","rfc.fit(X_train, y_train)\n","training_end = time.perf_counter()\n","prediction_start = time.perf_counter()\n","preds = rfc.predict(X_test)\n","prediction_end = time.perf_counter()\n","acc_rfc = (preds == y_test).sum().astype(float) / len(preds)*100\n","rfc_train_time = training_end-training_start\n","rfc_prediction_time = prediction_end-prediction_start\n","print(\"Scikit-Learn's Random Forest Classifier's prediction accuracy is: %3.2f\" % (acc_rfc))\n","print(\"Time consumed for training: %4.3f seconds\" % (rfc_train_time))\n","print(\"Time consumed for prediction: %6.5f seconds\" % (rfc_prediction_time))"]},{"cell_type":"markdown","metadata":{"_cell_guid":"15b03625-4960-45a9-9751-b8e3ea82cfed","_uuid":"ce81bb5c1480e8cf041153ffd04efa02fbf1c6ef"},"source":["#### Support Vector Machine Classifier"]},{"cell_type":"code","execution_count":null,"metadata":{"_cell_guid":"c50ce8dd-8a29-4cd5-80aa-68e13f2f6cea","_uuid":"546fec58bfc86bd223f59efeac49abc5041a91c0","collapsed":true,"scrolled":true},"outputs":[],"source":["svc = SVC()\n","training_start = time.perf_counter()\n","svc.fit(X_train, y_train)\n","training_end = time.perf_counter()\n","prediction_start = time.perf_counter()\n","preds = svc.predict(X_test)\n","prediction_end = time.perf_counter()\n","acc_svc = (preds == y_test).sum().astype(float) / len(preds)*100\n","svc_train_time = training_end-training_start\n","svc_prediction_time = prediction_end-prediction_start\n","print(\"Scikit-Learn's Support Vector Machine Classifier's prediction accuracy is: %3.2f\" % (acc_svc))\n","print(\"Time consumed for training: %4.3f seconds\" % (svc_train_time))\n","print(\"Time consumed for prediction: %6.5f seconds\" % (svc_prediction_time))"]},{"cell_type":"markdown","metadata":{"_cell_guid":"722e29dc-95a5-4a16-9120-50f43b15623c","_uuid":"7a87904d4510511e806c24db46c14b13dc24181d"},"source":["#### Dense Neural Network"]},{"cell_type":"code","execution_count":null,"metadata":{"_cell_guid":"bbe1b419-939f-42ed-b292-5239d2d20640","_uuid":"de3c554de2317463b3bdbc8fd0ebae1714c6bffc","collapsed":true,"scrolled":false},"outputs":[],"source":["# only output errors\n","tf.logging.set_verbosity(tf.logging.ERROR)\n","\n","# define initial hyperparameters\n","batch_size = 1\n","train_steps = 1000\n","hidden_units=[10,10]\n","\n","# create a numeric column for every column of the training set\n","feature_cols = []\n","\n","for key in X_train.keys():\n","     feature_cols.append(tf.feature_column.numeric_column(key=key))\n","\n","# build model with input function\n","dnn = tf.estimator.DNNClassifier(hidden_units=hidden_units, feature_columns=feature_cols, n_classes=3)\n","\n","\n","input_fn = tf.estimator.inputs.pandas_input_fn(x=X_train, \n","                                               y=y_train,\n","                                               batch_size=batch_size,\n","                                               shuffle=False)\n","\n","# fit the DNN to the training data\n","training_start = time.perf_counter()\n","dnn.train(input_fn = input_fn, steps=train_steps)\n","training_end = time.perf_counter()\n","\n","# predict test data\n","predict_input_func = tf.estimator.inputs.pandas_input_fn(x=X_test, shuffle=False)\n","\n","prediction_start = time.perf_counter()\n","preds = dnn.predict(input_fn=predict_input_func)\n","prediction_end = time.perf_counter()\n","\n","# evaluate results\n","final_preds = []\n","for pred in preds:\n","    final_preds.append(pred['class_ids'][0])\n","    \n","acc_dnn = (final_preds == y_test).sum().astype(float) / len(final_preds)*100\n","dnn_train_time = training_end-training_start\n","dnn_prediction_time = prediction_end-prediction_start\n","\n","print(\"Dense Neural Network Classifier's prediction accuracy is: %3.2f\" % (acc_dnn))\n","print(\"Time consumed for training: %4.3f seconds\" % (dnn_train_time))\n","print(\"Time consumed for prediction: %6.5f seconds\" % (dnn_prediction_time))"]},{"cell_type":"markdown","metadata":{"_cell_guid":"1d7e36a9-414c-4212-a86c-763ca8422b86","_uuid":"937ada3280f7e4de9b7cb022c47c0af1908238e6"},"source":["Let's compare the results. We will create a table for a more comprehensive overview."]},{"cell_type":"code","execution_count":null,"metadata":{"_cell_guid":"a79815ff-1d64-48ec-996d-32bcbba98250","_uuid":"c1fcf42a95967c7d25394c755e5014b027681d23","collapsed":true,"scrolled":false},"outputs":[],"source":["results = pd.DataFrame({\n","    'Model': ['KNN', 'Naive Bayes', \n","              'XGBoost', 'Random Forest', 'SVC', 'DNN'],\n","    'Score': [acc_knn, acc_gnb, acc_xgb, acc_rfc, acc_svc, acc_dnn],\n","    'Runtime Training': [knn_train_time, gnb_train_time, xgb_train_time, rfc_train_time, \n","                         svc_train_time, dnn_train_time],\n","    'Runtime Prediction': [knn_prediction_time, gnb_prediction_time, xgb_prediction_time, rfc_prediction_time,\n","                          svc_prediction_time, dnn_prediction_time]})\n","result_df = results.sort_values(by='Score', ascending=False)\n","result_df = result_df.set_index('Model')\n","result_df"]},{"cell_type":"markdown","metadata":{"_cell_guid":"07452e0e-4f68-4f7c-a368-97ed43251eb4","_uuid":"d108af400cc3ddb2bae13b870688d7149300c0f6"},"source":["We can see that both XGBoost and Scikit-Learn's Random Forest Classifier could achieve high accuracy.\n","\n","Gaussian Naive Bayes achieves less accuracy but needs a very little amount of time to both train and predict data.\n","\n","KNN performs about 10% worse than Naive Bayes.\n","\n","The Support Vector Machine Classifier has a bad accuracy, plus takes the a lot of time for its operations.\n","\n","The Dense Neural Network performs surprisingly unaccurate but takes very little time to predict."]},{"cell_type":"markdown","metadata":{"_cell_guid":"e1018b7a-65f5-46d3-9076-544fea2ea1ee","_uuid":"207493115d3d29e59fce048de6e25e82e7adc9fd"},"source":["### K Fold Cross Validation"]},{"cell_type":"markdown","metadata":{"_cell_guid":"58985896-8f98-4d6c-83e9-298b33a05057","_uuid":"a86a460ecd7d943d5f80d5b717030e306c6f6f18"},"source":["We will no perform k fold cross valdiation for the top 2 classifiers, i.e. XGBoost & Random Forest.\n","\n","We do this to get a more realistic result by testing the performance for 10 different train and test datasets and averaging the results. \n","\n","Cross validation ensures that the above result is not arbitary and gives a more reliable performance check."]},{"cell_type":"markdown","metadata":{"_cell_guid":"afaee1e4-5e36-405f-8bb6-82b3f90a8d8c","_uuid":"7c153977c26128dddc17304ad3e2a28c33f29c18"},"source":["#### Scikit-Learn's Random Forest Classifier"]},{"cell_type":"code","execution_count":null,"metadata":{"_cell_guid":"1c330298-5c44-4c4f-81c6-633082a644d0","_uuid":"277ceb9eea218726a798af0f9d1dfce6d56fbb0b","collapsed":true},"outputs":[],"source":["from sklearn.model_selection import cross_val_score\n","rfc_cv = RandomForestClassifier(n_estimators=100)\n","scores = cross_val_score(rfc_cv, X_train, y_train, cv=10, scoring = \"accuracy\")\n","print(\"Scores:\", scores)\n","print(\"Mean:\", scores.mean())\n","print(\"Standard Deviation:\", scores.std())"]},{"cell_type":"markdown","metadata":{"_cell_guid":"907bb075-8fb8-4627-806b-d6390e3bafea","_uuid":"cef974b7252639f4927f69a4ace837e9ee4a1811"},"source":["#### XGBoost"]},{"cell_type":"code","execution_count":null,"metadata":{"_cell_guid":"e12fade1-b4bb-4188-9e3a-8654acbd6627","_uuid":"3447d7d4313ce367d226fc2c111f9275eddb7307","collapsed":true,"scrolled":true},"outputs":[],"source":["xgb_cv = XGBClassifier(n_estimators=100)\n","scores = cross_val_score(xgb_cv, X_train, y_train, cv=10, scoring = \"accuracy\")\n","print(\"Scores:\", scores)\n","print(\"Mean:\", scores.mean())\n","print(\"Standard Deviation:\", scores.std())"]},{"cell_type":"markdown","metadata":{"_cell_guid":"80b3da72-6711-42fb-93e3-6670112cd96c","_uuid":"0941b7073b9d12b655476b6e44d453c065004478"},"source":["Cross validating the models showed that the accuracy values were in fact not arbitary and proofed that both models are performing very well. \n","\n","XGBoost showed a higher mean and lower standard deviation than the Scikit-Learn RFC.\n","\n","A high mean corresponds to a more stable performance and a low standard deviation corresponds to smaller range of results. "]},{"cell_type":"markdown","metadata":{"_cell_guid":"8d501e81-e7bf-4212-89b0-7e4dda3afd45","_uuid":"8538142325bbaf4c09449126856a8cf0725aeabe"},"source":["### Feature Importance"]},{"cell_type":"markdown","metadata":{"_cell_guid":"58abedc0-f5ba-497e-8037-4fd6e03e283a","_uuid":"d201ab385d30b19c444b79d0d0b7697f9161f888"},"source":["Decision Trees have the unique property of being able to order features by their ability to split between the classes.\n","\n","We will now visualize the features and their splitting ability."]},{"cell_type":"code","execution_count":null,"metadata":{"_cell_guid":"1684a815-20b7-47cd-ba41-d326249961e2","_uuid":"12d49c95b6b762d0cf43eba7c28818237fa3ae97","collapsed":true,"scrolled":true},"outputs":[],"source":["importances = pd.DataFrame({\n","    'Feature': X_train.columns,\n","    'Importance': xgb.feature_importances_\n","})\n","importances = importances.sort_values(by='Importance', ascending=False)\n","importances = importances.set_index('Feature')\n","importances"]},{"cell_type":"code","execution_count":null,"metadata":{"_cell_guid":"0852af00-7336-4b24-bc0a-2f2bb5f6454f","_uuid":"53754170a45188c7b5fcafc44c28707c850f8128","collapsed":true,"scrolled":false},"outputs":[],"source":["importances.plot.bar()"]},{"cell_type":"markdown","metadata":{"_cell_guid":"c20a7dc7-7226-4970-b480-44150db5ae8e","_uuid":"27931b2ee5325ff26b194ce646604116f315996b"},"source":["Here we can clearly see how PCA helped to improve the performance of our predictors as 2 of the principal components are in the top 3 features.\n","\n","The best (in terms of being able to split classes) is redshift.\n","\n","Mjd is the feature with the lowest importance during the classification process, we will therefore drop it from the dataframe."]},{"cell_type":"code","execution_count":null,"metadata":{"_cell_guid":"ea202e17-55a7-407f-b464-d618c31e64df","_uuid":"044570e7da58c81c04080f7d1290935741dd536b","collapsed":true},"outputs":[],"source":["sdss = sdss_df_fe.drop('mjd', axis=1)"]},{"cell_type":"markdown","metadata":{"_cell_guid":"c1ec976a-51d8-4bb2-b93b-bdf30981f277","_uuid":"87c2b7276cbe466417b8f696c93dc7055a3454ae"},"source":["### Summary"]},{"cell_type":"markdown","metadata":{"_cell_guid":"0d8dd596-537d-41be-81b2-df7c93aeddd5","_uuid":"e6ab1fa918694f2f4ca5cf93fa6c39550f4db3ec"},"source":["We trained different machine learning models to solve this classification problems. Without any further hyperparameter tuning XGBoost and Scikit-Learn's Random Forest Classifier performed the best.\n","\n","As XGBoost showed a little higher accuracy in most of the tests, we will continue only with this classifier."]},{"cell_type":"markdown","metadata":{"_cell_guid":"8a76e771-19a5-4bd0-ab08-9098f665d375","_uuid":"651c173faeaff4566ed8bff527731e56d3c8d684"},"source":["## XGBoost - Finding the best hyperparameters"]},{"cell_type":"markdown","metadata":{"_cell_guid":"1f618343-24e7-42f2-b95f-f7dac3502e2b","_uuid":"168f26a0042b88c070d287b2908c1b034e7b43e6"},"source":["Now it's time to look for the optimal hyperparameters - what does this mean?\n","\n","We will test our chosen model with different values for (almost) each of its tuning parameters and give back the parameters with which the model performed best.\n","\n","**The actual searching for optimal parameters is not done in this notebook since the operations can take some time and parallel editing would not be possible.\n","\n","We will write our transformed data set to disk so the tuning scripts can access it. "]},{"cell_type":"code","execution_count":null,"metadata":{"_cell_guid":"b95e6c43-7145-4735-8207-cb643359b8ee","_uuid":"d6ab85e0bd84d95871b2cd7e575541b1e1f9a3d4","collapsed":true},"outputs":[],"source":["sdss.to_csv('sdss_data.csv')"]},{"cell_type":"markdown","metadata":{"_cell_guid":"388d4fa4-69c8-47a4-b7c7-22863bdc544b","_uuid":"cebfc0611237c96beead9d52a5199412ece2451b"},"source":["The best parameters for prediction as found by the tuning tests are:\n","\n","- max_depth = 5\n","- min_child_weight = 1\n","- gamma = 0\n","- subsample = 0.8\n","- colsample_bytree = 0.8\n","- reg_alpha = 0.005"]},{"cell_type":"markdown","metadata":{"_cell_guid":"b77d76db-28e2-4f85-968c-bf1f81cf1117","_uuid":"97075d3a0b265ebb9dcc9b3a11aabecbcaf01222"},"source":["## XGBoost - Testing optimal hyperparameters"]},{"cell_type":"markdown","metadata":{"_cell_guid":"ed7d218a-b908-4e05-9842-3b800eb5e23a","_uuid":"64549d7e5f8df20030fc261b3cd10e86ef7a989a"},"source":["The following model implements the best performing model with optimal parameters evaluated by the hyperparameter tuning. We will expect the model to perform even better than before.\n","\n","Analytics Vidhya presented a really nice guide for tuning XGBoost. \n","\n","Please read more: https://www.analyticsvidhya.com/blog/2016/03/complete-guide-parameter-tuning-xgboost-with-codes-python/"]},{"cell_type":"code","execution_count":null,"metadata":{"_cell_guid":"32212a9d-b781-4028-af0b-6d473edd78eb","_uuid":"5014e6c1606f1773a70396c04714a32fa8641c94","collapsed":true},"outputs":[],"source":["X_train, X_test, y_train, y_test = train_test_split(sdss.drop('class', axis=1), sdss['class'],\n","                                                   test_size=0.33)"]},{"cell_type":"code","execution_count":null,"metadata":{"_cell_guid":"9a82daab-70d3-4641-8bdd-05f539bcac00","_uuid":"b62105d36b02c5bf23f357cdc642bad9ac7b5238","collapsed":true},"outputs":[],"source":["xgboost = XGBClassifier(max_depth=5, learning_rate=0.01, n_estimators=100, gamma=0, \n","                        min_child_weight=1, subsample=0.8, colsample_bytree=0.8, reg_alpha=0.005)\n","\n","xgboost.fit(X_train, y_train)\n","preds = xgboost.predict(X_test)\n","\n","accuracy = (preds == y_test).sum().astype(float) / len(preds)*100\n","\n","print(\"XGBoost's prediction accuracy WITH optimal hyperparameters is: %3.2f\" % (accuracy))"]},{"cell_type":"markdown","metadata":{"_cell_guid":"f96aeb1f-042f-4a68-9954-6f4b135f13fa","_uuid":"d87abb59478a4f75cf60ab9ac4f470ae04c78bae"},"source":["The parameter tuning did not improve the accuracy as excpected. We will therefore do a cross validation to test to get a more reliable result."]},{"cell_type":"code","execution_count":null,"metadata":{"_cell_guid":"b6caf968-61bd-4002-a596-48d59b2c1e82","_uuid":"971694ed539775cb89764b56b5a6de91b31a9773","collapsed":true,"scrolled":true},"outputs":[],"source":["xgb_cv = XGBClassifier(n_estimators=100)\n","scores = cross_val_score(xgb_cv, X_train, y_train, cv=10, scoring = \"accuracy\")\n","print(\"Scores:\", scores)\n","print(\"Mean:\", scores.mean())\n","print(\"Standard Deviation:\", scores.std())"]},{"cell_type":"markdown","metadata":{"_cell_guid":"f82c3f6d-e534-4428-a2e0-18ead8902907","_uuid":"1e473836eb2e64c3bfa77b3c580938b445345959"},"source":["Depending on the run the cross validation results vary from a little lower and a little better than before. \n","\n","This indicates that the parameter tuning was not as effective as expected - this could mean that XGBoost was actually close to its maximum performance capability on this data set.\n","\n","As we still have a good performance we will now continue with further evaluation of the performance of our model!"]},{"cell_type":"markdown","metadata":{"_cell_guid":"ab40ef62-a88c-40f2-8de1-557c1d3f3538","_uuid":"4145a1588a58c647a89277fe422f7336fed8e663"},"source":["## XGBoost - Evaluation"]},{"cell_type":"markdown","metadata":{"_cell_guid":"335223e5-543a-4220-b800-d6701fac8a0b","_uuid":"be081626a369bb71e249d2e8c1e908cb93be5278"},"source":["### Confusion Matrix"]},{"cell_type":"code","execution_count":null,"metadata":{"_cell_guid":"7de533e5-3f44-4c8e-a355-8506c091b3cb","_uuid":"dc5e8f82f108c866325cb2fafd8e819448454bf0","collapsed":true},"outputs":[],"source":["unique, counts = np.unique(sdss['class'], return_counts=True)\n","dict(zip(unique, counts))"]},{"cell_type":"code","execution_count":null,"metadata":{"_cell_guid":"979ba73d-780d-45e9-8be2-ccd51fdf4cec","_uuid":"edfebe01145e108c72ed26ac933050be270e8dff","collapsed":true},"outputs":[],"source":["predictions = cross_val_predict(xgb, sdss.drop('class', axis=1), sdss['class'], cv=3)\n","confusion_matrix(sdss['class'], predictions)"]},{"cell_type":"markdown","metadata":{"_cell_guid":"6caa5c87-1592-4a81-944d-2c95064f8c7b","_uuid":"fe9df9fdd4b5d0bd707cfb57022207ea4a280c06"},"source":["The first row shows that out of 4998 stars, **4962 were classified correctly as stars**. 29 stars were classified incorrectly as galaxies and 7 stars were classified incorrectly as quasars.\n","\n","The second row shows out of 850 quasars **827 were classified correctly**. 22 qsos were classified incorrectly as stars and 1 quasar was classified as galaxy.\n","\n","The last row tells us that out of 4152 galaxies **4147 were classified correctly.**. 5 galaxies were classified incorrectly as star.\n","\n","In total: \n","\n","We have only 64 objects which were classified incorrectly. Most of the objects were recognized as what they are."]},{"cell_type":"markdown","metadata":{"_cell_guid":"fed94503-ff9b-4ed8-b0d5-068ab9997c2e","_uuid":"620588ac4dd13964fed52434128fe04c4ca41648"},"source":["### Precision & Recall"]},{"cell_type":"code","execution_count":null,"metadata":{"_cell_guid":"70fd1280-a89b-4d82-b6bd-659be36ee9a2","_uuid":"658950ca13aa7d48c09e16a18e51d3748343f1d2","collapsed":true},"outputs":[],"source":["print(\"Precision:\", precision_score(sdss['class'], predictions, average='micro'))\n","print(\"Recall:\",recall_score(sdss['class'], predictions, average='micro'))"]},{"cell_type":"markdown","metadata":{"_cell_guid":"5e20bc1d-d511-4e4b-977f-c1dc03bdaff2","_uuid":"1cd5284547adaf76c3b54e722216b4eb38daf83e"},"source":["Precision is the fraction of events where the algorithm classified an object of type **t** correctly out of all occurences of the algorithm classifying objects of type **t**.\n","\n","Recall is the fraction of events where the algorithm classified an object of type **t** correctly when the true type of that object was actually **t**.\n","\n","Precision in our case:\n","\n","For every class its calculated how many objects were classified as stars (or galaxies or quasars) in relation to the amount of correct star (or galaxies or quasars) predictions. The results are averaged --> 99.36%.\n","\n","Recall in our case:\n","\n","For every class its calculated how many objects were classified as stars (or galaxies or quasars) in relation to the total amount of predictions where the object actually was a star (or galaxy or quasar). The results are averaged --> 99.36%.\n","\n","In both cases our algorithm did a very good job. The highest precision or recall value a predictor can have is 1.0."]},{"cell_type":"markdown","metadata":{"_cell_guid":"ea2c2e72-8733-4fe2-82f9-d545e67aee3b","_uuid":"2bb3f471d03c443abc684301f9cc86cf01b6fb00"},"source":["### F1-Score"]},{"cell_type":"markdown","metadata":{"_cell_guid":"891698cd-b105-4879-8ee3-f06c86415ae0","_uuid":"b3f24007d8dfd05ec3f24c231ecb933ec72f2cca"},"source":["One can combine precision and recall into one score, which is called the F-score. The F-score is computed with the harmonic mean of precision and recall. Note that it assigns much more weight to low values. As a result of that, the classifier will only get a high F-score, if both recall and precision are high."]},{"cell_type":"code","execution_count":null,"metadata":{"_cell_guid":"470445a5-a8fa-4fe1-a888-ec0ca77bb61a","_uuid":"ad9c64fb50ed540f5731a8230c0f2eb5791fc7f8","collapsed":true},"outputs":[],"source":["print(\"F1-Score:\", f1_score(sdss['class'], predictions, average='micro'))"]},{"cell_type":"markdown","metadata":{"_cell_guid":"64449fab-79a7-4f3f-a909-b39670e5fd5a","_uuid":"a63a0334223d4b644215d3f807ef99c8097394a2"},"source":["As precision and recall have the same value the F1-Score has automatically the same value too. Again, we are very close to 1.0 which indicates strong performance."]},{"cell_type":"markdown","metadata":{"_cell_guid":"9d6084f2-7121-4ef6-8f62-ce434f94f4ea","_uuid":"91e2f9ced455f5a8f02a80cb5ff035ba3228f7f4"},"source":["## Summary"]},{"cell_type":"markdown","metadata":{"_cell_guid":"93e6bbb6-c774-4ab9-b2a8-f71ee48196ed","_uuid":"d1f96b80f16c196f063143c7afc88254474dd4ea"},"source":["In this notebook learned how to get data from the SDSS RD14, analyze the data (we learned some very interesting facts about our space along the way), how to build a machine learning model to predict for unseen data from this data set and how to improve its performance (even though there was only a slight improvent).\n","We used XGBoost for predicting and evaluated its result.\n","\n","This project was very interesting to work on as I'm also interested in space and astronomy."]},{"cell_type":"code","execution_count":null,"metadata":{"_cell_guid":"a3720240-7114-4ad8-a67f-fdd4ffcca3dd","_uuid":"08db6990894ef8c612eeb4b9f4658e18bcfb6979","collapsed":true},"outputs":[],"source":[]}],"metadata":{"kernelspec":{"display_name":"Python 3","language":"python","name":"python3"},"language_info":{"codemirror_mode":{"name":"ipython","version":3},"file_extension":".py","mimetype":"text/x-python","name":"python","nbconvert_exporter":"python","pygments_lexer":"ipython3","version":"3.6.5"}},"nbformat":4,"nbformat_minor":1}