1_Clustering_and_Regression_Analysis.ipynb

{"cells":[{"metadata":{"_uuid":"3cad42cc96a72e367978375ec63c32da606c91a6"},"cell_type":"markdown","source":"<h1 align=\"center\"> The Impact of Medicine to your Wallet </h1>\n<img src=\"https://www.toonpool.com/user/1631/files/expensive_health_care_898485.jpg\" width=600>\n\n<h3>Before Starting: </h3>\nAs in other projects, your upvotes really mean a lot to me because it tells me that Kagglers are interested in the work I am proving to you guys. So I will appreciate if you could upvote this kernel if you enjoy the work I do. Looking to share some insights with Kagglers in the comment section. Also, if updates take longer than usual it is because of work at school nevertheless, I'll try to bring more interesting updates with regards to this project. Hope you enjoy the analysis!\n\n<h3> Brief Introduction </h3>\nIn this project, my main aim is to show ways to go deep into the data story-telling even though the dataset is small. Also, I will work on a model that could give us an approximation as to what will be the charges of the patients. Nevertheless, we must go deeply into what factors influenced the charge of a specific patient. In order to do this we must look for patterns in our data analysis and gain extensive insight of what the data is telling us.  Lastly, we will go step by step to understand the story behind the patients in this dataset only through this way we could have a better understanding of what features will help our model have a closer accuracy to the true patient charge. \n\n\n<h3>A Note to my Fellow Kagglers: </h3>\nI am still in the process of updating especially in the cluster analysis and multi regression phase. I am aware of the collinearity issue which I will update once I am done implementing and explaining the different types of cluster analysis we could use to determine specific groups within our population. Enjoy! <br><br>\n\n<h3>Table of Contents: </h3>\nI. [Importing Libraries](#import)<br>\nII. [Distribution of Medical Charges](#distribution)<br>\n\nIII. Age Analysis <br>\na) [Turning Age into Categorical Variables](#age_cat)<br>\nb) [Is there a relationship between BMI and Age?](#bmi_age)<br>\nc) [Comparing Independent Variables with ANOVA](#ANOVA)<br>\nd)[Who got charged more on Average by Age](#charged_age)<br><br>\n\nIV. Weight Status <br>\na) [Turning BMI into Categorical Variables](#bmi_cat) <br>\nb) [Weight Status vs Charges](#weight_charges) <br>\nc) [Obesity and the Impact of Smoking to the Wallet](#obese_smoker)<br>\nd)[ Distribution of Charges (Obese Smoker vs Obese non-Smoker)](#obesevsnonobese_smokers)<br>\ne) [Separation in Charges between Obese Smokers vs Non-Obese Smokers](#separation)<br><br>\n\nV. Regional Analysis:<br>\na) [Building a Contingency Table](#building_contingency)<br>\nb) [Average Patient Charge by Region](#average_region) <br>\nc) [Average charge by region depending on weight condition](#charge_condition) <br><br>\n\nVI. Unsupervised Learning Algorithms:<br>\na) [Performing Clustering Manually](#manual_cluster)<br>\nb) [K-Means Clustering](#Kmeans)<br>\nc) [Hierarchical Clustering](#Hierarchical_clustering)<br><br>\n\n\nVII. Correlations <br>\na) [Correlations and Bivariate Analysis](#correlations)<br>\nb) [Preprocessing Data and using Pipelines](#preprocess)<br><br>\n\n\nVIII. Statistics <br>\na) [Linear Models and Terms](#linear_models)<br>\nb) [Dealing with Collinearity](#collinearity)<br>\nc) [The importance of adding Smoking status to our model](#smoking_status)<br>\nd) [Dealing with Skewness](#skewness)<br>\ne) Using our Model in a symmetric distribution (To be updated)\n\n"},{"metadata":{"_uuid":"ea43cbf5edca54c14b0097d7ff1506cd6cef2f0d"},"cell_type":"markdown","source":"<a id=\"import\"></a>"},{"metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"import numpy as np\nimport pandas as pd\n\n\n# Plotly Packages\nfrom plotly import tools\nimport plotly.plotly as py\nimport plotly.figure_factory as ff\nimport plotly.graph_objs as go\nfrom plotly.offline import download_plotlyjs, init_notebook_mode, plot, iplot\ninit_notebook_mode(connected=True)\n\n# Matplotlib and Seaborn\nimport matplotlib.pyplot as plt\nimport seaborn as sns\nfrom string import ascii_letters\n\n# Statistical Libraries\nfrom scipy.stats import norm\nfrom scipy.stats import skew\nfrom scipy.stats.stats import pearsonr\nfrom scipy import stats\n\n\n# Regression Modeling\nimport statsmodels.api as sm\nfrom statsmodels.sandbox.regression.predstd import wls_prediction_std\n\n\n# Other Libraries\nimport warnings\nwarnings.filterwarnings(\"ignore\")","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"c24db7fa2a282ff39644edd6610d17d0ccca49fe","_kg_hide-input":true},"cell_type":"code","source":"df = pd.read_csv(\"../input/insurance.csv\")\ndf.head()\n\n# Let's store the original dataframe in another variable.\noriginal_df = df.copy()","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"1dea7192ffde3be479255574347d635dbbeab49a"},"cell_type":"markdown","source":"<h3> Distribution of Medical Charges </h3> \n<a id=\"distribution\"></a>\n<ul>\n    <li><b>Types of Distributions: </b> We have a <b> right skewed distribution </b> in which most patients are being charged between $2000 - $12000.</li>\n    <li><b>Using Logarithms:  </b> Logarithms helps us have a <b>normal distribution </b> which could help us in a number of different ways such as outlier detection, implementation of statistical concepts based on the central limit theorem and for our predictive modell in the foreseen future. (We will come to this later in the kernel) </li>\n    </ul>"},{"metadata":{"trusted":true,"_uuid":"18141c201905560777376b3ff616585b40c65475","_kg_hide-input":true},"cell_type":"code","source":"# Determine the distribution of charge\ncharge_dist = df[\"charges\"].values\nlogcharge = np.log(df[\"charges\"])\n\n\n\ntrace0 = go.Histogram(\n    x=charge_dist,\n    histnorm='probability',\n    name=\"Charges Distribution\",\n    marker = dict(\n        color = '#FA5858',\n    )\n)\ntrace1 = go.Histogram(\n    x=logcharge,\n    histnorm='probability',\n    name=\"Charges Distribution using Log\",\n    marker = dict(\n        color = '#58FA82',\n    )\n)\n\nfig = tools.make_subplots(rows=2, cols=1,\n                          subplot_titles=('Charge Distribution','Log Charge Distribution'),\n                         print_grid=False)\n\n\n\nfig.append_trace(trace0, 1, 1)\nfig.append_trace(trace1, 2, 1)\n\n\nfig['layout'].update(showlegend=True, title='Charge Distribution', bargap=0.05)\niplot(fig, filename='custom-sized-subplot-with-subplot-titles')","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"ce236bdcbd169d81205959cfbccf10a8bbc36526"},"cell_type":"markdown","source":"### Age Analysis:\n<a id=\"age_cat\"></a>\n<h4>Turning Age into Categorical Variables: </h4>\n<ul> \n    <li><b>Young Adult: </b> from 18 - 35 </li>\n    <li><b>Senior Adult: </b> from 36 - 55 </li>\n    <li><b>Elder: </b> 56 or older </li>\n    <li><b> Share of each Category: </b> Young Adults (42.9%), Senior Adults (41%) and Elder (16.1%) </li>\n    </ul>"},{"metadata":{"trusted":true,"_uuid":"3761d4571cc41a856bca688edfae1f991958203d","_kg_hide-input":true},"cell_type":"code","source":"df['age_cat'] = np.nan\nlst = [df]\n\nfor col in lst:\n    col.loc[(col['age'] >= 18) & (col['age'] <= 35), 'age_cat'] = 'Young Adult'\n    col.loc[(col['age'] > 35) & (col['age'] <= 55), 'age_cat'] = 'Senior Adult'\n    col.loc[col['age'] > 55, 'age_cat'] = 'Elder'\n    \n    \nlabels = df[\"age_cat\"].unique().tolist()\namount = df[\"age_cat\"].value_counts().tolist()\n\ncolors = [\"#ff9999\", \"#b3d9ff\", \" #e6ffb3\"]\n\ntrace = go.Pie(labels=labels, values=amount,\n               hoverinfo='label+percent', textinfo='value', \n               textfont=dict(size=20),\n               marker=dict(colors=colors, \n                           line=dict(color='#000000', width=2)))\n\ndata = [trace]\nlayout = go.Layout(title=\"Amount by Age Category\")\n\nfig = go.Figure(data=data, layout=layout)\niplot(fig, filename='basic_pie_chart')","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"a1d7af164f3c1e1d3548ae9bf579848e37b4fb3d"},"cell_type":"markdown","source":"<h3> Is there a Relationship between BMI and Age</h3>\n<a id=\"bmi_age\"></a>\n<ul>\n    <li><b>BMI frequency: </b> Most of the BMI frequency is concentrated between 27 - 33. </li>\n    <li><b> Correlations </b> Age and charges have a correlation of 0.29 while bmi and charges have a correlation of 0.19 </li>\n    <li><b> Relationship betweem BMI and Age: </b> The correlation for these two variables is 0.10 which is not that great. Therefore, we can disregard that age has a huge influence on BMI. </li>\n    </ul>"},{"metadata":{"trusted":true,"_uuid":"b6edb6e1f285aefc4222f13d949999a4d8d5c2bf","_kg_hide-input":true},"cell_type":"code","source":"bmi = [df[\"bmi\"].values.tolist()]\ngroup_labels = ['Body Mass Index Distribution']\n\ncolors = ['#FA5858']\n\nfig = ff.create_distplot(bmi, group_labels, colors=colors)\n# Add title\nfig['layout'].update(title='Normal Distribution <br> Central Limit Theorem Condition')\n\niplot(fig, filename='Basic Distplot')","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"0bab9ea41fac06a7746bafe22790ae600241b7cc","_kg_hide-input":true},"cell_type":"code","source":"corr = df.corr()\n\nhm = go.Heatmap(\n    z=corr.values,\n    x=corr.index.values.tolist(),\n    y=corr.index.values.tolist()\n)\n\n\ndata = [hm]\nlayout = go.Layout(title=\"Correlation Heatmap\")\n\nfig = dict(data=data, layout=layout)\niplot(fig, filename='labelled-heatmap')","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"1c7be215a60ae55f95ff48164648183b0ad2cecd","_kg_hide-input":true},"cell_type":"code","source":"young_adults = df[\"bmi\"].loc[df[\"age_cat\"] == \"Young Adult\"].values\nsenior_adult = df[\"bmi\"].loc[df[\"age_cat\"] == \"Senior Adult\"].values\nelders = df[\"bmi\"].loc[df[\"age_cat\"] == \"Elder\"].values\n\ntrace0 = go.Box(\n    y=young_adults,\n    name = 'Young Adults',\n    boxmean= True,\n    marker = dict(\n        color = 'rgb(214, 12, 140)',\n    )\n)\ntrace1 = go.Box(\n    y=senior_adult,\n    name = 'Senior Adults',\n    boxmean= True,\n    marker = dict(\n        color = 'rgb(0, 128, 128)',\n    )\n)\n\ntrace2 = go.Box(\n    y=elders,\n    name = 'Elders',\n    boxmean= True,\n    marker = dict(\n        color = 'rgb(247, 186, 166)',\n    )\n)\n\n\n\n\ndata = [trace0, trace1, trace2]\n\nlayout = go.Layout(title=\"Body Mass Index <br> by Age Category\", xaxis=dict(title=\"Age Category\", titlefont=dict(size=16)),\n                  yaxis=dict(title=\"Body Mass Index\", titlefont=dict(size=16)))\n\nfig = go.Figure(data=data, layout=layout)\niplot(fig)","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"dd59bc833e62eaf37e715bb894d1b7367fff7a65"},"cell_type":"markdown","source":"\n"},{"metadata":{"_uuid":"7c0effd289fe435fc90c52c65e032cbd388a3f50"},"cell_type":"markdown","source":"<h4> Comparing Independent Categorical Variables (ANOVA) </h4>\n<a id=\"ANOVA\"></a>\n<ul>\n    <li> <b> P-value: </b>The p-value being higher than 0.05 tells us that we take the Null hypothesis, meaning that there is no a significant change between the three age categories when it comes to Body Mass Index.</li>\n    </ul>"},{"metadata":{"trusted":true,"_uuid":"777321da108a782727a50d1d7c686877afd5c88d","_kg_hide-input":true},"cell_type":"code","source":"import statsmodels.api as sm\nfrom statsmodels.formula.api import ols\n\n\nmoore_lm = ols(\"bmi ~ age_cat\", data=df).fit()\nprint(moore_lm.summary())","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"5eded1df647e10cf64edf99ac96f73e5d8a5f0f2","_kg_hide-input":true},"cell_type":"code","source":"import plotly.plotly as py\nimport plotly.graph_objs as go\n\nya_smoker = df[\"bmi\"].loc[(df[\"age_cat\"] == \"Young Adult\") & (df[\"smoker\"] == \"yes\")].values\nsa_smoker = df[\"bmi\"].loc[(df[\"age_cat\"] == \"Senior Adult\") & (df[\"smoker\"] == \"yes\")].values\ne_smoker = df[\"bmi\"].loc[(df[\"age_cat\"] == \"Elder\") & (df[\"smoker\"] == \"yes\")].values\n\n# Non-Smokers\nya_nonsmoker = df[\"bmi\"].loc[(df[\"age_cat\"] == \"Young Adult\") & (df[\"smoker\"] == \"no\")].values\nsa_nonsmoker = df[\"bmi\"].loc[(df[\"age_cat\"] == \"Senior Adult\") & (df[\"smoker\"] == \"no\")].values\ne_nonsmoker = df[\"bmi\"].loc[(df[\"age_cat\"] == \"Elder\") & (df[\"smoker\"] == \"no\")].values\n\nx_data = ['Young A. Smoker', 'Young A. Non-Smoker',\n          'Senior A. Smoker', 'Senior A. Non-Smoker',\n          'Elder Smoker', 'Elder Non-Smoker',]\n\ny0 = ya_smoker\ny1 = ya_nonsmoker\ny2 = sa_smoker\ny3 = sa_nonsmoker\ny4 = e_smoker\ny5 = e_nonsmoker\n\ny_data = [y0,y1,y2,y3,y4,y5]\n\ncolors = ['rgba(251, 43, 43, 0.5)', 'rgba(125, 251, 137, 0.5)', \n          'rgba(251, 43, 43, 0.5)', 'rgba(125, 251, 137, 0.5)', \n          'rgba(251, 43, 43, 0.5)', 'rgba(125, 251, 137, 0.5)']\n\ntraces = []\n\nfor xd, yd, cls in zip(x_data, y_data, colors):\n        traces.append(go.Box(\n            y=yd,\n            name=xd,\n            boxpoints='all',\n            jitter=0.5,\n            whiskerwidth=0.2,\n            fillcolor=cls,\n            marker=dict(\n                size=2,\n            ),\n            line=dict(width=1),\n        ))\n\nlayout = go.Layout(\n    title='Body Mass Index of Smokers Status by Age Category',\n    xaxis=dict(\n    title=\"Status\",\n    titlefont=dict(\n    size=16)),\n    yaxis=dict(\n        title=\"Body Mass Index\",\n        autorange=True,\n        showgrid=True,\n        zeroline=True,\n        dtick=5,\n        gridcolor='rgb(255, 255, 255)',\n        gridwidth=1,\n        zerolinecolor='rgb(255, 255, 255)',\n        zerolinewidth=2,\n        titlefont=dict(\n        size=16)\n    ),\n    margin=dict(\n        l=40,\n        r=30,\n        b=80,\n        t=100,\n    ),\n    paper_bgcolor='rgb(255, 255, 255)',\n    plot_bgcolor='rgb(255, 243, 192)',\n    showlegend=False\n)\n\nfig = go.Figure(data=traces, layout=layout)\niplot(fig)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"edf258b59f8ba1c2b807e6b727bfc8424afff644"},"cell_type":"markdown","source":"### Who got charged more on Average by Age?\n<a id=\"charged_age\"></a>\n<ul> \n    <li><b>Patient Charge Mean: </b> For <b>young adults </b> it is 7,944, for <b> Senior Adults </b> it is 14,785  and for the <b>elder</b> it is 18,795. </li>\n     <li><b>Patient Charge Median: </b> For <b>young adults </b> it is 4,252, for <b> Senior Adults </b> it is 9,565  and for the <b>elder</b> it is 13,429. </li>\n    <li> <b> Mean and the Median:</b> Sometimes we must be careful when using the mean since it is prone to be affected by outliers.</li>\n    </ul>"},{"metadata":{"trusted":true,"_uuid":"536c613a2a1c381fe5bbebc902d0b00b318a5585","_kg_hide-input":true},"cell_type":"code","source":"# Mean could be affected easily by outliers or extreme cases.\n# Means\navg_ya_charge = df[\"charges\"].loc[df[\"age_cat\"] == \"Young Adult\"].mean()\navg_sa_charge = df[\"charges\"].loc[df[\"age_cat\"] == \"Senior Adult\"].mean()\navg_e_charge = df[\"charges\"].loc[df[\"age_cat\"] == \"Elder\"].mean()\n\n# Median\nmed_ya_charge = df[\"charges\"].loc[df[\"age_cat\"] == \"Young Adult\"].median()\nmed_sa_charge = df[\"charges\"].loc[df[\"age_cat\"] == \"Senior Adult\"].median()\nmed_e_charge = df[\"charges\"].loc[df[\"age_cat\"] == \"Elder\"].median()\n\naverage_plot = go.Bar(\n    x=['Young Adults', 'Senior Adults', 'Elder'],\n    y=[avg_ya_charge, avg_sa_charge, avg_e_charge],\n    name='Mean',\n    marker=dict(\n        color=\"#F5B041\"\n    )\n)\nmed_plot = go.Bar(\n    x=['Young Adults', 'Senior Adults', 'Elder'],\n    y=[med_ya_charge, med_sa_charge, med_e_charge],\n    name='Median',\n    marker=dict(\n        color=\"#48C9B0\"\n    )\n)\n\n\nfig = tools.make_subplots(rows=1, cols=2, specs=[[{}, {}]],\n                          subplot_titles=('Average Charge by Age','Median Charge by Age'),\n                         shared_yaxes=True, print_grid=False)\n\n\nfig.append_trace(average_plot, 1, 1)\nfig.append_trace(med_plot, 1, 2)\n\n\nfig['layout'].update(showlegend=True, title='Age Charges', xaxis=dict(title=\"Age Category\"), yaxis=dict(title=\"Patient Charges\"), bargap=0.15)\niplot(fig, filename='custom-sized-subplot-with-subplot-titles')","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"6e098c3cc4bb6dd0e025da1ac3dce9a1fc5072e4"},"cell_type":"markdown","source":"### Weight Status:\nhttps://www.cancer.org/cancer/cancer-causes/diet-physical-activity/body-weight-and-cancer-risk/adult-bmi.html\n\n<h4>Turning BMI into Categorical Variables: </h4>\n<a id=\"bmi_cat\"></a>\n<ul>\n    <li><b>Under Weight: </b> Body Mass Index (BMI) $<$ 18.5  </li>\n    <li><b>Normal Weight: </b>  Body Mass Index (BMI) $\\geq$ 18.5 and  Body Mass Index (BMI) $<$ 24.9 </li>\n    <li><b>Overweight: </b>  Body Mass Index (BMI) $\\geq$ 25 and  Body Mass Index (BMI) $<$ 29.9</li>\n    <li><b>Obese: </b>  Body Mass Index (BMI) $>$ 30 </li>\n    </ul>"},{"metadata":{"trusted":true,"_uuid":"0100ae1f9fcc2294860b0d8c29f876d2c18d832f","_kg_hide-input":true},"cell_type":"code","source":"df[\"weight_condition\"] = np.nan\nlst = [df]\n\nfor col in lst:\n    col.loc[col[\"bmi\"] < 18.5, \"weight_condition\"] = \"Underweight\"\n    col.loc[(col[\"bmi\"] >= 18.5) & (col[\"bmi\"] < 24.986), \"weight_condition\"] = \"Normal Weight\"\n    col.loc[(col[\"bmi\"] >= 25) & (col[\"bmi\"] < 29.926), \"weight_condition\"] = \"Overweight\"\n    col.loc[col[\"bmi\"] >= 30, \"weight_condition\"] = \"Obese\"\n    \ndf.head()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"4bd0deba0dd09ac627ea6194b6b863a6ec879966","_kg_hide-input":true},"cell_type":"code","source":"# Create subpplots\nf, (ax1, ax2, ax3) = plt.subplots(ncols=3, figsize=(18,8))\n\n# I wonder if the cluster that is on the top is from obese people\nsns.stripplot(x=\"age_cat\", y=\"charges\", data=df, ax=ax1, linewidth=1, palette=\"Reds\")\nax1.set_title(\"Relationship between Charges and Age\")\n\n\nsns.stripplot(x=\"age_cat\", y=\"charges\", hue=\"weight_condition\", data=df, ax=ax2, linewidth=1, palette=\"Set2\")\nax2.set_title(\"Relationship of Weight Condition, Age and Charges\")\n\nsns.stripplot(x=\"smoker\", y=\"charges\", hue=\"weight_condition\", data=df, ax=ax3, linewidth=1, palette=\"Set2\")\nax3.legend_.remove()\nax3.set_title(\"Relationship between Smokers and Charges\")\n\nplt.show()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"6d58d34154b464a8cd6455cd639103098ab98226"},"cell_type":"code","source":"# Make sure we don't have any null values\ndf[df.isnull().any(axis=1)]","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"b5bfd92c34651e6ea02cb52f93259d0015dadc00"},"cell_type":"markdown","source":"<h3> Weight Status vs Charges</h3>\n<a id=\"weight_charges\"></a>\n\n<ul>\n    <li> <b>Overweight: </b> Notice how there are two groups of people that get significantly charged more than the other group of overweight people. </li>\n    <li><b>Obese: </b> Same thing goes with the obese group, were a significant group is charged more than the other group. </li>\n    </ul>\n    \n**    Hmmm, might the smoking status have to do with this phenomenon?**<br> Let's find out!"},{"metadata":{"trusted":true,"_uuid":"5b04bcd46dcb732f8c11e49ba62388eeec4008ae","_kg_hide-input":true},"cell_type":"code","source":"fig = ff.create_facet_grid(\n    df,\n    x='age',\n    y='charges',\n    color_name='weight_condition',\n    show_boxes=False,\n    marker={'size': 10, 'opacity': 1.0},\n    colormap={'Underweight': 'rgb(208, 246, 130)', 'Normal Weight': 'rgb(166, 246, 130)',\n             'Overweight': 'rgb(251, 232, 238)', 'Obese': 'rgb(253, 45, 28)'}\n)\n251, 232, 238\n\n\nfig['layout'].update(title=\"Weight Status vs Charges\", width=800, height=600, plot_bgcolor='rgb(251, 251, 251)', \n                     paper_bgcolor='rgb(255, 255, 255)')\n\n\niplot(fig, filename='facet - custom colormap')","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"f4ef78746793d14889479cc4d958a479aee1e232","_kg_hide-input":true},"cell_type":"code","source":"# First find the average or median of the charges obese people paid.\n\nobese_avg = df[\"charges\"].loc[df[\"weight_condition\"] == \"Obese\"].mean()\n\ndf[\"charge_status\"] = np.nan\nlst = [df]\n\n\nfor col in lst:\n    col.loc[col[\"charges\"] > obese_avg, \"charge_status\"] = \"Above Average\"\n    col.loc[col[\"charges\"] < obese_avg, \"charge_status\"] = \"Below Average\"\n    \ndf.head()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"24e090f4c1d2edcbac5452ece0efe49b235558c8","_kg_hide-input":true},"cell_type":"code","source":"# No Nulls\ndf[\"charge_status\"].isnull().sum()","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"9ba9b92023caf0a4315bf4eba9212d429d9a907f"},"cell_type":"markdown","source":"<h3> Obesity and the Impact of Smoking to the Wallet: </h3>\n<a id=\"obese_smoker\"> </a>\n<ul> \n    <li> Notice in the <b> charges </b> box how smoking looks to have a certain impact on medical costs.</li>\n    <li> Let's find out how much of a difference there is between the group of <b> obese </b> patients that smoke compared to the group of <b> obese </b> patients that don't smoke.</li>\n    </ul>"},{"metadata":{"trusted":true,"_uuid":"a01f702d12a684401f0c4de0a714d9ac054caf6c","_kg_hide-input":true},"cell_type":"code","source":"import seaborn as sns\nsns.set(style=\"ticks\")\npal = [\"#FA5858\", \"#58D3F7\"]\n\nsns.pairplot(df, hue=\"smoker\", palette=pal)\nplt.title(\"Smokers\")","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"83c5ce81dcb4fc96c48cb202aa201e908609b870","_kg_hide-input":true},"cell_type":"code","source":"# What Percentage of Obese that Smoked Paid aBove Average from the total obese patients?\n# 79% of Obese were non-smokers while the 21% left were smokers\ntotal_obese = len(df.loc[df[\"weight_condition\"] == \"Obese\"])\n\nobese_smoker_prop = len(df.loc[(df[\"weight_condition\"] == \"Obese\") & (df[\"smoker\"] == \"yes\")])/total_obese\nobese_smoker_prop = round(obese_smoker_prop, 2)\n\nobese_nonsmoker_prop = len(df.loc[(df[\"weight_condition\"] == \"Obese\") & (df[\"smoker\"] == \"no\")])/total_obese\nobese_nonsmoker_prop = round(obese_nonsmoker_prop, 2)\n\n\n# Average charge by obese_smokers and obese_nonsmoker\ncharge_obese_smoker = df.loc[(df[\"weight_condition\"] == \"Obese\") & (df[\"smoker\"] == \"yes\")].mean()\ncharge_obese_nonsmoker = df.loc[(df[\"weight_condition\"] == \"Obese\") & (df[\"smoker\"] == \"no\")].mean()","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"1f3a5e40aea004444659779785cad0712fa1f1d4"},"cell_type":"markdown","source":"<h4>Distribution of Charges (Obese Smoker vs Obese non-Smoker) </h4>\n<a id=\"obesevsnonobese_smokers\"></a>\n<ul> \n    <li><b>Violin Plots: </b> We will be using violin plots to <b>compare the distributions of patients of the obese group who are smokers and non-smokers. </b> </li>\n    <li><b> Obese smokers distribution: </b> Most obese smokers pay around 40k in medical costs! </li>\n    <li> <b>Obese non-smokers distribution: </b> Most obese non-smokers pay 8k in medical costs. </li>\n    <li> <b> Smoking a factor to the wallet </b> Smoking is defninitely a big factor for obese patients when it comes to medical cost. A difference of more than 30k!</li>\n    </ul>"},{"metadata":{"trusted":true,"_uuid":"9e2af84d5685c0358c55d67fe7c73383e375943f","_kg_hide-input":true},"cell_type":"code","source":"\npointspossmoker = [-0.9,-1.1,-0.6,-0.3]\npointposnonsmoker = [0.45,0.55,1,0.4]\nshowLegend = [True,False,False,False]\n\ndata = []\nfor i in range(0,len(pd.unique(df['weight_condition']))):\n    male = {\n            \"type\": 'violin',\n            \"x\": df['weight_condition'][ (df['smoker'] == 'yes') & (df['weight_condition'] == pd.unique(df['weight_condition'])[i]) ],\n            \"y\": df['charges'][ (df['smoker'] == 'yes') & (df['weight_condition'] == pd.unique(df['weight_condition'])[i]) ],\n            \"legendgroup\": 'Smoker',\n            \"scalegroup\": 'Smoker',\n            \"name\": 'Smoker',\n            \"side\": 'negative',\n            \"box\": {\n                \"visible\": True\n            },\n            \"points\": 'all',\n            \"pointpos\": pointspossmoker[i],\n            \"jitter\": 0,\n            \"scalemode\": 'count',\n            \"meanline\": {\n                \"visible\": True\n            },\n            \"line\": {\n                \"color\": '#DF0101'\n            },\n            \"marker\": {\n                \"line\": {\n                    \"width\": 2,\n                    \"color\": '#F78181'\n                }\n            },\n            \"span\": [\n                0\n            ],\n            \"showlegend\": showLegend[i]\n        }\n    data.append(male)\n    female = {\n            \"type\": 'violin',\n            \"x\": df['weight_condition'] [ (df['smoker'] == 'no') & (df['weight_condition'] == pd.unique(df['weight_condition'])[i]) ],\n            \"y\": df['charges'] [ (df['smoker'] == 'no') & (df['weight_condition'] == pd.unique(df['weight_condition'])[i]) ],\n            \"legendgroup\": 'Non-Smoker',\n            \"scalegroup\": 'Non-Smoker',\n            \"name\": 'Non-Smoker',\n            \"side\": 'positive',\n            \"box\": {\n                \"visible\": True\n            },\n            \"points\": 'all',\n            \"pointpos\": pointposnonsmoker[i],\n            \"jitter\": 0,\n            \"scalemode\": 'count',\n            \"meanline\": {\n                \"visible\": True\n            },\n            \"line\": {\n                \"color\": '#00FF40'\n            },\n            \"marker\": {\n                \"line\": {\n                    \"width\": 2,\n                    \"color\": '#81F781'\n                }\n            },\n            \"span\": [\n                0\n            ],\n            \"showlegend\": showLegend[i]\n        }\n    data.append(female)\n        \n\nfig = {\n    \"data\": data,\n    \"layout\" : {\n        \"title\": \"Charges Distribution of Obese Patients<br><i>Group by Smoking Status\",\n        \"yaxis\": {\n            \"zeroline\": False,\n            \"title\": \"Patient Charges\",\n            \"titlefont\": {\n                \"size\": 16\n            }\n        },\n        \"violingap\": 0,\n        \"violingroupgap\": 0,\n        \"violinmode\": \"overlay\"\n    }\n}\n\n\niplot(fig, filename='violin/advanced', validate = False)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"0540e3fca71cc134cffff0a6b8be6d2dec08a7be","_kg_hide-input":true},"cell_type":"code","source":"# Hmmm we have to look closer into Obsese there is an obvious difference\n\nchargedist_sm = df[\"charges\"].loc[(df[\"weight_condition\"] == \"Obese\") & (df[\"smoker\"] == \"yes\")].values\nchargedist_nsm = df[\"charges\"].loc[(df[\"weight_condition\"] == \"Obese\") & (df[\"smoker\"] == \"no\")].values\n\ntrace0 = go.Box(\n    y=chargedist_sm,\n    name = 'Obese Smokers',\n    marker = dict(\n        color = '#DF0101',\n    )\n)\ntrace1 = go.Box(\n    y=chargedist_nsm,\n    name = 'Obese Non-Smokers',\n    marker = dict(\n        color = '#00FF40',\n    )\n)\n\n\ndata = [trace0, trace1]\n\n\nlayout = dict(title=\"Deeper Look into Obese condition by Smoking status\",\n             xaxis=dict(\n             title=\"Status\",\n             titlefont=dict(\n             size=16)),\n             yaxis=dict(title=\"Patient Charges\",\n                       titlefont=dict(size=16)),\n              )\n\nfig = go.Figure(data=data, layout=layout)\niplot(fig)","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"e8bbca189ea741d6a51ab4b221588faa65143e4b"},"cell_type":"markdown","source":"<h4>Separation in Charges between Obese Smokers vs Non-Obese Smokers </h4>\n<a id=\"separation\"></a>\nIn this chart we can visualize how can <b> separate</b> obese smokers and obese non-smokers into different <b> clusters </b> of groups. Therefore, we can say that smoking is a <b>characteristic</b> that definitely affects patient's charges."},{"metadata":{"trusted":true,"_uuid":"e2cfc01e41283c9485ad64b2096a16f8d2f5b928","_kg_hide-input":true},"cell_type":"code","source":"# Create a Scatter Plot with all the Obese\nobese_smoker = df.loc[(df[\"weight_condition\"] == \"Obese\") & (df[\"smoker\"] == \"yes\")]\nobese_nonsmoker = df.loc[(df[\"weight_condition\"] == \"Obese\") & (df[\"smoker\"] == \"no\")]\n\n\ntrace0 = go.Scatter(\n    x = obese_smoker[\"age\"].values,\n    y = obese_smoker[\"charges\"].values,\n    name = 'Smokers',\n    mode = 'markers',\n    marker = dict(\n        size = 10,\n        color = '#DF0101',\n        line = dict(\n            width = 2,\n            color = 'rgb(0, 0, 0)'\n        )\n    )\n)\n\ntrace1 = go.Scatter(\n    x = obese_nonsmoker[\"age\"].values,\n    y = obese_nonsmoker[\"charges\"].values,\n    name = 'Non-Smokers',\n    mode = 'markers',\n    marker = dict(\n        size = 10,\n        color = '#00FF40',\n        line = dict(\n            width = 2,\n        )\n    )\n)\n\ndata = [trace0, trace1]\n\nlayout = dict(title = 'Clear Separation between Obese Smokers and Non-Smokers in Charges',\n              yaxis = dict(zeroline = False,\n                          title=\"Patient Charges\",\n                          titlefont=dict(size=16)),\n              xaxis = dict(zeroline = False,\n                          title=\"Age of the Patient\",\n                          titlefont=dict(\n                          size=16))\n             )\n\nfig = dict(data=data, layout=layout)\niplot(fig, filename='styled-scatter')","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"ce5e6db5c5b2ea741159590197de516aeb25130c"},"cell_type":"markdown","source":"<h2 align=\"center\"><b>Regional Analysis: </b></h2>\n<h3> Building a Contingency Table: </h3>\n<a id=\"building_contingency\">"},{"metadata":{"trusted":true,"_uuid":"58fb21d954701a70cc37e8038bb30f24d1b5b301","_kg_hide-input":true},"cell_type":"code","source":"# Cont table by weight condition and sex\ndf.head()\n\nregion_smoker = pd.crosstab(df['smoker'], df['region']).apply(lambda x: x/x.sum() * 100)\nregion_smoker","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"5a4c091b1517caecaa4c1278d91e0a712afcdcb3","_kg_hide-input":true},"cell_type":"code","source":"plt.style.use('seaborn-whitegrid')\n\nticks = df['age_cat'].unique()\ncolors = ['#ff2424', '#90ee90']\n\n\nax = sns.catplot(x=\"age_cat\", y=\"charges\", hue=\"smoker\",\n                 col=\"region\", aspect=.6,\n                 kind=\"swarm\", palette=colors, data=df);\n\nax.set_xticklabels(labels = ticks, rotation=45)\n\nplt.show()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"c8974c0df0884339df18fbbe771f03f6869805d1","_kg_hide-input":true},"cell_type":"code","source":"# Obesity per region\nregion_weight = pd.crosstab(df['weight_condition'], df['region']).apply(lambda x: x/x.sum() * 100)\nregion_weight = round(region_weight, 2)\nregion_weight","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"4fbec3e4392eed86e9d93130c9cbd94098b04020"},"cell_type":"markdown","source":"<h3>Average Patient Charge by Region: </h3>\n<a id=\"average_region\"></a>\n<ul>\n    <li><b> Median Patient Charges: </b> The <b>NorthEast</b> is the region that pays the most on average while the <b>SouthWest</b> is the one that pays </li>\n    <li><b>Obese group: </b> From the obese group, the Southwest is the region where obese patients pay the most..</li>\n    <li><b>Overweight: </b> From the obese group, the NorthWest is the region where obese patients pay the most.  </li>\n    <li><b>Normal Weight: </b> From the obese group, the SouthEast is the region where obese patients pay the most.  </li>\n    <li><b>Underweight: </b> From the obese group, the NorthWest is the region where obese patients pay the most.</li>\n    </ul>"},{"metadata":{"trusted":true,"_uuid":"11288a683829602dbd821002d7650fc14d727520","_kg_hide-input":true},"cell_type":"code","source":"df.head()\n\n\n# Average charge by Region\ndf[\"region\"].unique()\n\n# Median Charges per Region\nsouthwest = np.median(df[\"charges\"].loc[df[\"region\"] == \"southwest\"].values)\nsoutheast = np.median(df[\"charges\"].loc[df[\"region\"] == \"southeast\"].values)\nnorthwest = np.median(df[\"charges\"].loc[df[\"region\"] == \"northwest\"].values)\nnortheast = np.median(df[\"charges\"].loc[df[\"region\"] == \"northeast\"].values)\n\nlst = [southwest, southeast, northwest, northeast]\n\ndata = [go.Scatterpolar(\n  r = [southwest, southeast, northwest, northeast],\n  theta = ['SouthWest', 'SouthEast', 'NorthWest', 'NorthEast'],\n  fill = 'toself'\n)]\n\nlayout = go.Layout(\n    title=\"Median Charged to Patients by Region\",\n    paper_bgcolor = \"rgb(255, 255, 224)\",\n  polar = dict(\n    radialaxis = dict(\n      visible = True,\n      range = [0, max(lst)]\n    )\n  ),\n  showlegend = False\n)\n\n\nfig = go.Figure(data=data, layout=layout)\niplot(fig, filename = \"radar/basic\")","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"05cf2312b514e37157254621fa73d75c7843eb95"},"cell_type":"markdown","source":"<h4>Average Charge by Region depending on the Weight Condition: </h4>\n<a id=\"charge_condition\"></a>\n"},{"metadata":{"trusted":true,"_uuid":"9cd67393099c207ed75993b8a72cd3866ad89587","_kg_hide-input":true},"cell_type":"code","source":"# Weight Condition by Region Radar plots\n\ndf[\"weight_condition\"].unique()\n\n# Average charges for overweight patients by region \nsw_overweight = np.mean(df[\"charges\"].loc[(df[\"region\"] == \"southwest\") & (df[\"weight_condition\"] == \"Overweight\")].values)\nse_overweight = np.mean(df[\"charges\"].loc[(df[\"region\"] == \"southeast\") & (df[\"weight_condition\"] == \"Overweight\")].values)\nnw_overweight = np.mean(df[\"charges\"].loc[(df[\"region\"] == \"northwest\") & (df[\"weight_condition\"] == \"Overweight\")].values)\nne_overweight = np.mean(df[\"charges\"].loc[(df[\"region\"] == \"northeast\") & (df[\"weight_condition\"] == \"Overweight\")].values)\n\n# Obese\nsw_obese = np.mean(df[\"charges\"].loc[(df[\"region\"] == \"southwest\") & (df[\"weight_condition\"] == \"Obese\")].values)\nse_obese = np.mean(df[\"charges\"].loc[(df[\"region\"] == \"southeast\") & (df[\"weight_condition\"] == \"Obese\")].values)\nnw_obese = np.mean(df[\"charges\"].loc[(df[\"region\"] == \"northwest\") & (df[\"weight_condition\"] == \"Obese\")].values)\nne_obese = np.mean(df[\"charges\"].loc[(df[\"region\"] == \"northeast\") & (df[\"weight_condition\"] == \"Obese\")].values)\n\n# Normal Weight\nsw_nw = np.mean(df[\"charges\"].loc[(df[\"region\"] == \"southwest\") & (df[\"weight_condition\"] == \"Normal Weight\")].values)\nse_nw = np.mean(df[\"charges\"].loc[(df[\"region\"] == \"southeast\") & (df[\"weight_condition\"] == \"Normal Weight\")].values)\nnw_nw = np.mean(df[\"charges\"].loc[(df[\"region\"] == \"northwest\") & (df[\"weight_condition\"] == \"Normal Weight\")].values)\nne_nw = np.mean(df[\"charges\"].loc[(df[\"region\"] == \"northeast\") & (df[\"weight_condition\"] == \"Normal Weight\")].values)\n\n# Underweight\nsw_uw = np.mean(df[\"charges\"].loc[(df[\"region\"] == \"southwest\") & (df[\"weight_condition\"] == \"Underweight\")].values)\nse_uw = np.mean(df[\"charges\"].loc[(df[\"region\"] == \"southeast\") & (df[\"weight_condition\"] == \"Underweight\")].values)\nnw_uw = np.mean(df[\"charges\"].loc[(df[\"region\"] == \"northwest\") & (df[\"weight_condition\"] == \"Underweight\")].values)\nne_uw = np.mean(df[\"charges\"].loc[(df[\"region\"] == \"northeast\") & (df[\"weight_condition\"] == \"Underweight\")].values)\n\n# Labels\nweight_labels = df[\"weight_condition\"].unique().tolist()\n\n# List per weight condition\nsw_weights = [sw_overweight, sw_obese, sw_nw, sw_uw]\nse_weights = [se_overweight, se_overweight, se_nw, se_uw]\nnw_weights = [nw_overweight, nw_overweight, nw_nw, nw_uw]\nne_weights = [ne_overweight, ne_overweight, ne_nw, ne_uw]\n\ndata = [\n    go.Scatterpolar(\n        mode=\"lines+markers\",\n        r = sw_weights,\n        theta = weight_labels,\n        fill = 'toself',\n        name=\"SouthWest\",\n        line=dict(\n            color=\"rgba(0, 128, 128, 0.95)\"\n        ),\n        marker=dict(\n            color=\"rgba(0, 74, 147, 1)\",\n            symbol=\"square\",\n            size=8\n        ),\n        subplot = \"polar\"\n    ),\n    go.Scatterpolar(\n        mode=\"lines+markers\",\n        r = se_weights,\n        theta = weight_labels,\n        fill = 'toself',\n        name=\"SouthEast\",\n        line=dict(\n            color=\"rgba(255, 72, 72, 0.95)\"\n        ),\n        marker=dict(\n            color=\"rgba(219, 0, 0, 1)\",\n            symbol=\"square\",\n            size=8\n        ),\n        subplot = \"polar2\"\n    ),\n    go.Scatterpolar(\n        mode=\"lines+markers\",\n        r = nw_weights,\n        theta = weight_labels,\n        fill = 'toself',\n        name=\"NorthWest\",\n        line=dict(\n            color=\"rgba(72, 255, 72, 0.95)\"\n        ),\n        marker=dict(\n            color=\"rgba(0, 147, 74, 1)\",\n            symbol=\"square\",\n            size=8\n        ),\n        subplot = \"polar3\"\n    ),\n       go.Scatterpolar(\n        mode=\"lines+markers\",\n        r = ne_weights,\n        theta = weight_labels,\n        fill = 'toself',\n        name=\"NorthEast\",\n        line=dict(\n            color=\"rgba(247, 133, 11, 0.95)\"\n        ),\n        marker=dict(\n            color=\"rgba(245, 168, 86, 1)\",\n            symbol=\"square\",\n            size=8\n        ),\n        subplot = \"polar4\"\n    )\n]\n\nlayout = go.Layout(\n    title=\"Average Patient Charges <br> by Region <br>(Depending on the Patient's Weight Condition)\",\n    showlegend = False,\n     paper_bgcolor = \"rgb(252, 234, 161)\",\n    polar = dict(\n      domain = dict(\n        x = [0, 0.46],\n        y = [0.56, 1]\n      ),\n      radialaxis = dict(\n        tickfont = dict(\n          size = 6\n        )\n      ),\n      angularaxis = dict(\n        tickfont = dict(\n          size = 8\n        ),\n        rotation = 40,\n        direction = \"clockwise\"\n      )\n    ),\n    polar2 = dict(\n      domain = dict(\n        x = [0, 0.46],\n        y = [0, 0.44]\n      ),\n      radialaxis = dict(\n        tickfont = dict(\n          size = 6\n        )\n      ),\n      angularaxis = dict(\n        tickfont = dict(\n          size = 8\n        ),\n        rotation = 40,\n        direction = \"clockwise\"\n      ),\n    ),\n    polar3 = dict(\n      domain = dict(\n       x = [0.54, 1],\n        y = [0.56, 1]\n      ),\n      radialaxis = dict(\n        tickfont = dict(\n          size = 6\n        )\n      ),\n      angularaxis = dict(\n        tickfont = dict(\n          size = 8\n        ),\n        rotation = 40,\n        direction = \"clockwise\"\n      ),\n    ),\n        polar4 = dict(\n      domain = dict(\n        x = [0.54, 1],\n        y = [0, 0.44]\n      ),\n      radialaxis = dict(\n        tickfont = dict(\n          size = 6\n        )\n      ),\n      angularaxis = dict(\n        tickfont = dict(\n          size = 8\n        ),\n        rotation = 40,\n        direction = \"clockwise\"\n      ),\n    ))\n\nfig = go.Figure(data=data, layout=layout)\nfig['layout'].update(height=800, width=800)\niplot(fig, filename='polar/directions')","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"b04773901ddd54994769cf7bbe53a78efa303b04"},"cell_type":"markdown","source":"<h2 align=\"center\">Unsupervised Learning: </h2>\n<h3>Performing Clustering in a Manual Way: </h3>\n<a id=\"manual_cluster\"></a>\nIn the first plot we will do a cluster analysis in a manual form and see what our eyes can discover. Here are the following results from the manual cluster analysis performed.\n<ul>\n    <li><b> Age and Charges: </b> We can see there is a slight increase in charges depending on the age of the patient. </li>\n    <li><b> Obese Clusters: </b> We can see that for each age group there are clusters of the obese group in the top part of charges. </li>\n    <li><b>Are these clusters Smokers? </b> As seen in the right chart, most of this clusters are definitely smokers.  </li>\n    \n</ul>\n\nThis is somewhat the same as we have explored in the previous section, but I wanted to add the factor of <b> age </b> to see until what extent this variable patient's charges. Nevertheless, we can still confirm that being obese and a smoker is the <b>grand major </b> factor at least for the group of obese and overweight patients. <br><br>\n\n<b> Thanks to  Bhumika Bhatt for adding the factor of age with the combination of weight condition and smoking status </b>"},{"metadata":{"trusted":true,"_uuid":"467a4e2b53ad6ac53a5f1f900d01077af192cc79","_kg_hide-input":true},"cell_type":"code","source":"# Two subplots one with weight condition and the other with smoker.\n\nf, (ax1, ax2) = plt.subplots(ncols=2, figsize=(18,8))\nsns.scatterplot(x=\"bmi\", y=\"charges\", hue=\"weight_condition\", data=df, palette=\"Set1\", ax=ax1)\nax1.set_title(\"Relationship between Charges and BMI by Weight Condition\")\nax1.annotate('Obese Cluster \\n (Does this cluster has \\n the Smoking Attribute?)', xy=(37, 50000), xytext=(30, 60000),\n            arrowprops=dict(facecolor='black'),\n            fontsize=12)\nsns.scatterplot(x=\"bmi\", y=\"charges\", hue=\"smoker\", data=df, palette=\"Set1\", ax=ax2)\nax2.set_title(\"Relationship between Charges and BMI by Smoking Condition\")\nax2.annotate('Obese Smoker Cluster ', xy=(35, 48000), xytext=(20, 60000),\n            arrowprops=dict(facecolor='black'),\n            fontsize=12)\nax2.annotate('The Impact of Smoking to \\n Charges on other \\n Weight Conditions ', xy=(25, 26000), xytext=(17, 40000),\n            arrowprops=dict(facecolor='black'),\n            fontsize=12)","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"dedef806ecb6799b4971611016ffdbb6b0cd5ba0"},"cell_type":"markdown","source":"<h2 align=\"center\"> <b>Understanding KMeans Clustering: </b> </h2>\n<a id=\"Kmeans\"></a>\n<img src=\"https://media.giphy.com/media/12vVAGkaqHUqCQ/giphy.gif\"><br><br>\n\n\n<h3> Before We Start Explaining the Elbow Method and KMeans Neighbors.  </h3>\nThe elbow method is mostly used in unsupervised learning algorithms to determine the <b> optimal number of clusters </b> that should be used to find specific unknown groups within our population. We used the <b>yellowbrick </b> library to implement a simple elbow method and to determine the appropiate number of clusters in our KMeans algorithm. <br><br>\n\n<h3>Terms to Know: </h3>\n<ul>\n    <li><b> Cluster Centroids: </b> The cluster centroid is the most representative point of a specific cluster. So, if we decide to find three clusters, we will have three cluster centroid. </li>\n    <li><b>Euclidean Distance: </b> Is the distance between two data points and this term is essential when gathering the distance between the cluster centroids and the data points.  </li>\n    <li><b> Elbow Method: </b> The elbow method is a technique used to choose the most optimal number of clusters. Remember, in Kmeans clustering we add the number of clusters in a manual way, so the elbow method is useful when using Kmeans. Why is it called the Elbow method? Because as more iterations run to find the optimal number of clusters, the line will take the shape of the arm and the optimal number of clusters is the point that is in the elbow part of the arm. </li>\n</ul>\n\n<h3> Euclidean Distance Formula: </h3>\n$\\LARGE d(x,y) = \\sqrt{(x_1 - x_2)^2 + (y_1 - y_2)^2}$<br>\n\nwhere: <br>\n$x_1$ = X-axis value of  data observation <br>\n$x_2$ = X-axis value of the cluster centroid <br>\n$y_1$ = Y-axis value of the data observation. <br>\n$y_2$ = Y-axis value of the cluster centroid. <br>\n\n<h3> Explaining the Elbow Method: </h3>\nThe elbow method finds the average sum of squares distance between the cluster centroid and the data observations. As the number of cluster increases the average sum of squares decreases. Basically, as the number of clusters increases, the distance between the data points and the centroids decreases as well. Whenever, we see the \"elbow\" that is a rule of thumb to consider the optimal number of clusters.\n\n\n<h3 align=\"center\">How the Elbow Method Looks Like: </h3>\n<img src=\"https://media.licdn.com/dms/image/C4E12AQHnQ7zAvCnvZg/article-inline_image-shrink_1500_2232/0?e=1553126400&v=beta&t=_QjaAH2ydyzjy3qVwAvV9Zi1WictiLOEam9L6YuSdI4\" width=600>\n\n"},{"metadata":{"trusted":true,"_uuid":"7aed64f0a6be551c149c74c6d0acb0f90eca8dc2","_kg_hide-input":true},"cell_type":"code","source":"from sklearn.cluster import KMeans\nfrom yellowbrick.cluster import KElbowVisualizer\n\nfig = plt.figure(figsize=(12,8))\n\n# KNears Neighbors \ndf.head()\noriginal_df.head()\n\nX = df[[\"bmi\", \"charges\"]]\n\n\n# Instantiate the clustering model and visualizer\nmodel = KMeans()\nvisualizer = KElbowVisualizer(model, k=(2,6))\n\nvisualizer.fit(X)    # Fit the data to the visualizer\nvisualizer.poof()  ","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"7301104e13d3109dada8479f6077451466b8275f","_kg_hide-input":true},"cell_type":"code","source":"import pandas as pd\nfrom sklearn.datasets import load_iris\nfrom sklearn.cluster import KMeans\nimport matplotlib.pyplot as plt\n\nkmeans = KMeans(n_clusters=3)  \nkmeans.fit(X)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"fec2ddbf93883a4b974e67a152dc8a50ef3a5bf0","_kg_hide-input":true},"cell_type":"code","source":"# Printing the Centroids\nprint(kmeans.cluster_centers_)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"6c85a46f6bcb7e4d832cccd0dcb4cf40234a32b8","_kg_hide-input":true},"cell_type":"code","source":"print(kmeans.labels_)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"7386b98eb15e85d65a6709ce447b3479d9029786","_kg_hide-input":true},"cell_type":"code","source":"fig = plt.figure(figsize=(12,8))\n\nplt.scatter(X.values[:,0], X.values[:,1], c=kmeans.labels_, cmap=\"Set1_r\", s=25)\nplt.scatter(kmeans.cluster_centers_[:,0] ,kmeans.cluster_centers_[:,1], color='black', marker=\"x\", s=250)\nplt.title(\"Kmeans Clustering \\n Finding Unknown Groups in the Population\", fontsize=16)\nplt.show()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"f7794edf613351c5aa6d8b9420723c26a1062a6b"},"cell_type":"markdown","source":"<h2> Adding more Clusters with Hierarchical Clustering:</h2>\n<img src=\"https://www.saedsayad.com/images/Clustering_h1.png\" width=800><br><br>\n\nWe use a dendogram as show above to show how two clusters are merged into one big cluster.\n\n<h3> Two types of Approaches: </b>\n<ul>\n    <li><b> Agglomerative (bottom-up):</b> Each observation starts as <b>one cluster.</b> Based on the distance of those clusters (in this case observations) it will <b> merge </b> into one cluster. For instance, let's say observation A (Cluster A) and observation B (Cluster B) are two different clusters that are within a close distance. So in essence, it combines the two nearest clusters into one bigger cluster. Remember, as with K-Means clustering each cluster is represented by the <b>centroid </b> which is the average position of the data points (observations). Then each centroid will merge either with other centroids from other clusters or with individual observations which are considered to be individual clusters.</li>\n    <li><b>Divisive (Top-Bottom): </b>  With this approach we start at the top with <b>one big cluster.</b> The cluster will be partitioned at a point where it splits the big cluster into two big ones and it will run K-means into each of the clusters splitting the data further down. This will get to a point where the observations cannot be split any more since each observation becomes its own cluster. If you want to understand more about K-means look at the previos example at the top of hierarchical clustering. In practice, the divisive method is not used as often as the Agglomerative method. </li>\n    \n</ul>\n    \n<h3> <b>What does the dendogram represent?</b> </h3>\nTo make things simple starting from the bottom, each leaf means that two observations or clusters have been merged into one bigger cluster like we explained in the Agglomerative approach. In other words it shows the \"hierarchical\" relationship between the clusters.\n"},{"metadata":{"trusted":true,"_uuid":"987d1b23b07fe4b7cddbc1414240098b923ba39b","_kg_hide-input":true},"cell_type":"code","source":"from sklearn.cluster import AgglomerativeClustering\n\nX = df[[\"bmi\", \"charges\"]]\n\nagglomerative_clustering = AgglomerativeClustering(n_clusters=4).fit(X)\nagglomerative_clustering","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"d98eeb3eecc27e4973fa93853aa30bd85d2ab34c","_kg_hide-input":true},"cell_type":"code","source":"from scipy.cluster.hierarchy import dendrogram, linkage\n\n# 5% of the data \nsample_df = df.sample(frac=.05)\n\nsample_X = sample_df[[\"bmi\", \"charges\"]]\n\nsample_agglomerative_clustering = AgglomerativeClustering(n_clusters=4).fit(sample_X)\nsample_agglomerative_clustering\n\n\nlinked = linkage(sample_agglomerative_clustering.children_, 'single')","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"dd537fef7c553b4126ffacb88274f2e2d25a59cf","_kg_hide-input":true},"cell_type":"code","source":"agglomerative_clustering.labels_","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"b62906f6599027deabd67918e5cd2e2bd0728aba","_kg_hide-input":true},"cell_type":"code","source":"print(plt.style.available)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"63dc486ce52eb6d637437d3d33d3d955f8e2f4a7","_kg_hide-input":true},"cell_type":"code","source":"plt.style.use(\"Solarize_Light2\")\n\nfig, (ax1, ax2) = plt.subplots(ncols=2, figsize=(16,6))\n\nax1.scatter(X.values[:,0], X.values[:,1], c=agglomerative_clustering.labels_, cmap=\"Set1_r\", s=25)\nax1.set_title(\"Agglomerative Clustering\", fontsize=16)\n\ndendrogram(linked,  \n            orientation='top',\n            labels=sample_agglomerative_clustering.labels_,\n            distance_sort='descending',\n            show_leaf_counts=False,\n          ax=ax2)\n\nax2.set_title(\"Dendogram on Agglomerative Clustering\")\n\nplt.show()","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"1ca84c59142226da5f9968781973bb59f910e122"},"cell_type":"markdown","source":"<h3> "},{"metadata":{"_uuid":"480615f5ffdcfc1ef22c2924898f3ee2259d1beb"},"cell_type":"markdown","source":"<h3> Stress Levels (Will elaborate further on this new feature) </h3>\n"},{"metadata":{"trusted":true,"_uuid":"8a04ef5e9cfdb8aefae3f0803a48f0ad9973de77","_kg_hide-input":true},"cell_type":"code","source":"# BMI + Children / Age\ndf[\"stress_level\"] = df[\"children\"] * df[\"age\"] / df[\"bmi\"]\n\ndf.head()","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"83fa4e1fac3d584b668d30fe96923fe289b4c5ce"},"cell_type":"markdown","source":"<h3> Correlations and Bivariate Analysis: </h3>\n<a id=\"correlations\"></a>"},{"metadata":{"trusted":true,"_uuid":"1724ba7389cb307460b015118532536fca2de385","_kg_hide-input":true},"cell_type":"code","source":"corr = df.corr()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"d23af4877bd14f9f4662f13698996c11b44780d8","_kg_hide-input":true},"cell_type":"code","source":"\n\nmask = np.zeros_like(corr, dtype=np.bool)\nmask[np.triu_indices_from(mask)] = True\n\n# Set up the matplotlib figure\nf, ax = plt.subplots(figsize=(12, 10))\n\n# Generate a custom diverging colormap\ncmap = sns.diverging_palette(220, 10, as_cmap=True)\n\n# Draw the heatmap with the mask and correct aspect ratio\nsns.heatmap(corr, mask=mask, cmap=cmap, vmax=.3, center=0,\n            square=True, linewidths=.5, cbar_kws={\"shrink\": .5})\n\nplt.title(\"Diagonal Correlation Matrix\", fontsize=20)\n\nplt.show()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"06136935732a6c9697f8eb32457a5c4586619987","_kg_hide-input":true},"cell_type":"code","source":"# pd.scatter_matrix(df, figsize=(12, 8))\n\n\nsns.pairplot(df,hue=\"charge_status\")\nplt.show()","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"b0edaae9c71c53d8b3f1865fe4d6db25c429f0be"},"cell_type":"markdown","source":"<h3> Preprocessing our Data and using Pipelines: </h3>\n<a id=\"preprocess\"></a>\nIn this section we will preprocess our data using what we call <b> Pipelines </b>. So why do we use pipelines? Although, in this example pipelines might not be necessarry if you want to automate various processes for instance, scaling the features, replacing missing values with the median and other processes pipelines are fantastic because they automatize the whole process.We will also use the class CategoricalEncoder from Scikit-Learn, they should already have it as a class but in this case I just copied the class from Scikit-Learn's library. Let's see how it works!"},{"metadata":{"trusted":true,"_uuid":"af00b9e6c2970839b8ddc7c2210fa47c8a47a297","_kg_hide-input":false,"_kg_hide-output":false},"cell_type":"code","source":"# In this section we will preprocess our data\n# First we should split our original data.\n\nfrom sklearn.model_selection import train_test_split\n\n# Shuffle our dataset before splitting\n\noriginal_df = original_df.sample(frac=1, random_state=1)\n\nX = original_df.drop(\"charges\", axis=1)\ny = original_df[\"charges\"]\n\n# Split into both training and testing\nX_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.2)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"9d512dec968400c38e67dd91b78003e0522551db","_kg_hide-input":false},"cell_type":"code","source":"from sklearn.preprocessing import StandardScaler\nfrom sklearn.pipeline import Pipeline\n\nfrom sklearn.base import BaseEstimator, TransformerMixin\nfrom sklearn.utils import check_array\nfrom sklearn.preprocessing import LabelEncoder\nfrom scipy import sparse\n\nclass CategoricalEncoder(BaseEstimator, TransformerMixin):\n    \"\"\"Encode categorical features as a numeric array.\n    The input to this transformer should be a matrix of integers or strings,\n    denoting the values taken on by categorical (discrete) features.\n    The features can be encoded using a one-hot aka one-of-K scheme\n    (``encoding='onehot'``, the default) or converted to ordinal integers\n    (``encoding='ordinal'``).\n    This encoding is needed for feeding categorical data to many scikit-learn\n    estimators, notably linear models and SVMs with the standard kernels.\n    Read more in the :ref:`User Guide <preprocessing_categorical_features>`.\n    Parameters\n    ----------\n    encoding : str, 'onehot', 'onehot-dense' or 'ordinal'\n        The type of encoding to use (default is 'onehot'):\n        - 'onehot': encode the features using a one-hot aka one-of-K scheme\n          (or also called 'dummy' encoding). This creates a binary column for\n          each category and returns a sparse matrix.\n        - 'onehot-dense': the same as 'onehot' but returns a dense array\n          instead of a sparse matrix.\n        - 'ordinal': encode the features as ordinal integers. This results in\n          a single column of integers (0 to n_categories - 1) per feature.\n    categories : 'auto' or a list of lists/arrays of values.\n        Categories (unique values) per feature:\n        - 'auto' : Determine categories automatically from the training data.\n        - list : ``categories[i]`` holds the categories expected in the ith\n          column. The passed categories are sorted before encoding the data\n          (used categories can be found in the ``categories_`` attribute).\n    dtype : number type, default np.float64\n        Desired dtype of output.\n    handle_unknown : 'error' (default) or 'ignore'\n        Whether to raise an error or ignore if a unknown categorical feature is\n        present during transform (default is to raise). When this is parameter\n        is set to 'ignore' and an unknown category is encountered during\n        transform, the resulting one-hot encoded columns for this feature\n        will be all zeros.\n        Ignoring unknown categories is not supported for\n        ``encoding='ordinal'``.\n    Attributes\n    ----------\n    categories_ : list of arrays\n        The categories of each feature determined during fitting. When\n        categories were specified manually, this holds the sorted categories\n        (in order corresponding with output of `transform`).\n    Examples\n    --------\n    Given a dataset with three features and two samples, we let the encoder\n    find the maximum value per feature and transform the data to a binary\n    one-hot encoding.\n    >>> from sklearn.preprocessing import CategoricalEncoder\n    >>> enc = CategoricalEncoder(handle_unknown='ignore')\n    >>> enc.fit([[0, 0, 3], [1, 1, 0], [0, 2, 1], [1, 0, 2]])\n    ... # doctest: +ELLIPSIS\n    CategoricalEncoder(categories='auto', dtype=<... 'numpy.float64'>,\n              encoding='onehot', handle_unknown='ignore')\n    >>> enc.transform([[0, 1, 1], [1, 0, 4]]).toarray()\n    array([[ 1.,  0.,  0.,  1.,  0.,  0.,  1.,  0.,  0.],\n           [ 0.,  1.,  1.,  0.,  0.,  0.,  0.,  0.,  0.]])\n    See also\n    --------\n    sklearn.preprocessing.OneHotEncoder : performs a one-hot encoding of\n      integer ordinal features. The ``OneHotEncoder assumes`` that input\n      features take on values in the range ``[0, max(feature)]`` instead of\n      using the unique values.\n    sklearn.feature_extraction.DictVectorizer : performs a one-hot encoding of\n      dictionary items (also handles string-valued features).\n    sklearn.feature_extraction.FeatureHasher : performs an approximate one-hot\n      encoding of dictionary items or strings.\n    \"\"\"\n\n    def __init__(self, encoding='onehot', categories='auto', dtype=np.float64,\n                 handle_unknown='error'):\n        self.encoding = encoding\n        self.categories = categories\n        self.dtype = dtype\n        self.handle_unknown = handle_unknown\n\n    def fit(self, X, y=None):\n        \"\"\"Fit the CategoricalEncoder to X.\n        Parameters\n        ----------\n        X : array-like, shape [n_samples, n_feature]\n            The data to determine the categories of each feature.\n        Returns\n        -------\n        self\n        \"\"\"\n\n        if self.encoding not in ['onehot', 'onehot-dense', 'ordinal']:\n            template = (\"encoding should be either 'onehot', 'onehot-dense' \"\n                        \"or 'ordinal', got %s\")\n            raise ValueError(template % self.handle_unknown)\n\n        if self.handle_unknown not in ['error', 'ignore']:\n            template = (\"handle_unknown should be either 'error' or \"\n                        \"'ignore', got %s\")\n            raise ValueError(template % self.handle_unknown)\n\n        if self.encoding == 'ordinal' and self.handle_unknown == 'ignore':\n            raise ValueError(\"handle_unknown='ignore' is not supported for\"\n                             \" encoding='ordinal'\")\n\n        X = check_array(X, dtype=np.object, accept_sparse='csc', copy=True)\n        n_samples, n_features = X.shape\n\n        self._label_encoders_ = [LabelEncoder() for _ in range(n_features)]\n\n        for i in range(n_features):\n            le = self._label_encoders_[i]\n            Xi = X[:, i]\n            if self.categories == 'auto':\n                le.fit(Xi)\n            else:\n                valid_mask = np.in1d(Xi, self.categories[i])\n                if not np.all(valid_mask):\n                    if self.handle_unknown == 'error':\n                        diff = np.unique(Xi[~valid_mask])\n                        msg = (\"Found unknown categories {0} in column {1}\"\n                               \" during fit\".format(diff, i))\n                        raise ValueError(msg)\n                le.classes_ = np.array(np.sort(self.categories[i]))\n\n        self.categories_ = [le.classes_ for le in self._label_encoders_]\n\n        return self\n\n    def transform(self, X):\n        \"\"\"Transform X using one-hot encoding.\n        Parameters\n        ----------\n        X : array-like, shape [n_samples, n_features]\n            The data to encode.\n        Returns\n        -------\n        X_out : sparse matrix or a 2-d array\n            Transformed input.\n        \"\"\"\n        X = check_array(X, accept_sparse='csc', dtype=np.object, copy=True)\n        n_samples, n_features = X.shape\n        X_int = np.zeros_like(X, dtype=np.int)\n        X_mask = np.ones_like(X, dtype=np.bool)\n\n        for i in range(n_features):\n            valid_mask = np.in1d(X[:, i], self.categories_[i])\n\n            if not np.all(valid_mask):\n                if self.handle_unknown == 'error':\n                    diff = np.unique(X[~valid_mask, i])\n                    msg = (\"Found unknown categories {0} in column {1}\"\n                           \" during transform\".format(diff, i))\n                    raise ValueError(msg)\n                else:\n                    # Set the problematic rows to an acceptable value and\n                    # continue `The rows are marked `X_mask` and will be\n                    # removed later.\n                    X_mask[:, i] = valid_mask\n                    X[:, i][~valid_mask] = self.categories_[i][0]\n            X_int[:, i] = self._label_encoders_[i].transform(X[:, i])\n\n        if self.encoding == 'ordinal':\n            return X_int.astype(self.dtype, copy=False)\n\n        mask = X_mask.ravel()\n        n_values = [cats.shape[0] for cats in self.categories_]\n        n_values = np.array([0] + n_values)\n        indices = np.cumsum(n_values)\n\n        column_indices = (X_int + indices[:-1]).ravel()[mask]\n        row_indices = np.repeat(np.arange(n_samples, dtype=np.int32),\n                                n_features)[mask]\n        data = np.ones(n_samples * n_features)[mask]\n\n        out = sparse.csc_matrix((data, (row_indices, column_indices)),\n                                shape=(n_samples, indices[-1]),\n                                dtype=self.dtype).tocsr()\n        if self.encoding == 'onehot-dense':\n            return out.toarray()\n        else:\n            return out","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"fe0373d7d146b5682b66e0ceb2ae151ed1cdecad","_kg_hide-input":false,"_kg_hide-output":false},"cell_type":"code","source":"from sklearn.base import BaseEstimator, TransformerMixin\n\n# A class to select numerical or categorical columns \n# since Scikit-Learn doesn't handle DataFrames yet\nclass DataFrameSelector(BaseEstimator, TransformerMixin):\n    def __init__(self, attribute_names):\n        self.attribute_names = attribute_names\n    def fit(self, X, y=None):\n        return self\n    def transform(self, X):\n        return X[self.attribute_names]","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"e010932e6198724f58562d92e8bfbc9e1b7bc402"},"cell_type":"code","source":"from sklearn.pipeline import Pipeline\nfrom sklearn.preprocessing import StandardScaler\nfrom sklearn.pipeline import FeatureUnion\n\n# Children as categorical (ordinal varibale)\nX_train[\"children\"] = X_train[\"children\"].astype(\"object\") \n\n# Separate numerics and categorical values\nnumerics = X_train.select_dtypes(exclude=\"object\")\ncategoricals = X_train.select_dtypes(include=\"object\")\n\n# Pipelines\nnumerical_pipeline = Pipeline([\n    (\"select_numeric\", DataFrameSelector(numerics.columns.tolist())),\n    (\"std_scaler\", StandardScaler()),\n])\n\ncategorical_pipeline =  Pipeline([\n    (\"select_numeric\", DataFrameSelector(categoricals.columns.tolist())),\n    (\"std_scaler\", CategoricalEncoder(encoding=\"onehot-dense\")),\n])\n\nmain_pipeline = FeatureUnion(transformer_list=[\n    ('num_pipeline', numerical_pipeline),\n    ('cat_pipeline', categorical_pipeline)\n])\n\n# Scale our features from our training data\nscaled_xtrain = main_pipeline.fit_transform(X_train)","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"cb42e088598c65c501e5e6aa09bf7117e2ee4863"},"cell_type":"markdown","source":"<h2 align=\"center\">  Statistics: </h2>\n<h3>Linear Models: </h3>\n<a id=\"linear_models\"></a>\nIn this section I used the statsmodels library to share with the community some important statistical measures I will give a brief definition to the concepts in which I consider are the most important in this case. <br>\n<b> Terms: </b>\n<ul> \n    <li><b>Adjusted R-Squared: </b> This tells us how close are our data points to the regression line (our predictions). The higher the adjusted squared, the higher the accuracy of our prediction line. </li>\n    <li><b>Skewness:</b> This tells us whether the distribution is symetric (at 0 or closer to 0), Right-Skewed which means that most of the values are to the left (Skewness statistic becomes larger) and Left skewed most observations are to the right (this means that the skewness statistic is in the negative zone.) </li>\n    <li><b> Kurtosis: </b> Is a statistical measure that helps us understand how skew is the distribution and the peakedness of our data. In this case we have a kurtosis greater than zero which means that the peak is wider with heavier tails, meaning there are some outliers in the observation. </li>\n    <li> <b>Skew:</b>  Our skewness is greater than zero which means that it is positive skew. Meaning that most observations are concentrated in the low section of the distribution.  </li>\n    <li><b>Collinearity: </b> Collinearity basically means that when two independent variables are highly correlated, it is hard to determine the real coefficients on each feature. Why is this important? Well, if we want to know how much weight that independent variable has towards predicting a specific value colinearity should be avoided.   </li>\n    </ul>\n    \n    \n  <h3>More on R-Square and Adjusted R-Square: </h3>  \n  <ul>\n    <li><b> What does R-Square assumes?</b> That every independent variable in the model explains the different levels of variations in the dependent variable (including the redundant independent variables.)  The more variables you add to the model, the more the R-Square will increase independently of how much significance the variable adds to the model.</li>\n    <li><b>Why Adjusted-R Squared?</b> To make things simple, adjusted-R Squared only uses the variables that adds significance to our model. This is one of the main reasons why we used adjusted R-squared rather than R-squared as a measure of performance of our model.  </li>\n</ul>\n    \n<h3>More on Collinearity: </h3>\n<ul> \n    <li><b>Does Collinearity affect the final value of our predictions?</b> The answer is no (well at least most of the time does not affect the value of the predictions) however,  collinearity prevents us from understanding the true coefficients of the independent variables which is important if we want to better understand how much weight the independent variables have towards a dependent variable. </li>\n    <li><b> What do we do to avoid collinearity?</b> First step, is to determine which variables are highly correlated and then decide which independent variable to choose for inclusion into our regression model. </li>\n    <li><b>Why is determining collinearity important?</b> When we find a high correlation between two features, this assumes that they are measuring the same label (output). So removing one of those features does not have a negative impact in the performance of our model and it makes our model simpler. Collinearity is solved by dropping those features that are highly correlated. </li>\n  </ul> \n  \n  \n <h3>Formulas: </h3>\n<b> R-Squared: </b><br>\n$\\LARGE R^2 = 1 - \\frac{SS Regression} {SS Total}$ <br><br>\nwhere: <br>\nSS Regression = Sum Squared Regression Error <br>\nSS Total = Sum Squared Total Error  <br><br>\n\n<b>Adjusted-R Squared: </b><br>\n$\\LARGE R^2 adjusted = 1 - \\frac{(1 - R^2) - (N - 1)}{N - p - 1}$ <br><br>\n\nwhere: <br>\n$R^2$ = Sample R-Squared <br>\n$p$ = Number of Predictors <br>\n$N$ = Total Sample Size <br><br>\n\n<b> Skeweness: </b><br>\nSkewness = $\\LARGE\\frac{\\sum(x_i - \\overline{x})}{(n - 1)v^\\frac{3}{2}}$<br><br>\n\nwhere: <br>\n$\\large V = \\frac{\\sum(x_i - \\overline{x})^2}{(n - 1)}$<br>\n$\\large x $= Predictor variable <br>\n$\\large n $= Number of Values<br>\n$\\large\\overline {x} $= Sample Mean $\n\n\n    \n <b>A Special Thanks Note: </b> A special thanks to <b> Piyush Rastogi </b>for bringing the collinearity issue to my attention. "},{"metadata":{"trusted":true,"_uuid":"5b69cae45fc7776d02579e36762a8a2790f256eb","_kg_hide-input":false},"cell_type":"code","source":"# Let's create the training set by combining the previous X_train and y_train.\ntrain = X_train.join(y_train, lsuffix='_X_train', rsuffix='_y_train')\ntest = X_test.join(y_test, lsuffix='_X_test', rsuffix='_y_test')\n\n# Random seed\nnp.random.seed(42)\n\n# Shuffle Randomly the training set\ntrain = train.sample(frac=1)\ntrain.head()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"2ecfc3fb314122594a33821ae5180491ee88cbc8","_kg_hide-input":true},"cell_type":"code","source":"X_train = sm.add_constant(scaled_xtrain)\ny_train = y_train.values\n\nmodel = sm.OLS(y_train, X_train)\nresults = model.fit()\nprint(results.summary())","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"5e748463e0160bee429ad07d4a78b6a4ed70121f","_kg_hide-input":true},"cell_type":"code","source":"print('Parameters: ', results.params)\nprint('R2: ', results.rsquared)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"eb64c0051e1ee6acd640fead5a155fb301040144","_kg_hide-input":true},"cell_type":"code","source":"res = sm.OLS(y_train, X_train).fit()\nprint(res.summary())","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"7e3b16fdfd3918bdb47c278edb15556ff5be5759","_kg_hide-input":true},"cell_type":"code","source":"print('Parameters: ', res.params)\nprint('Standard errors: ', res.bse)\nprint('Predicted values: ', res.predict())","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"814155a40f94795db45a7c9206cd239ac23ee1ea"},"cell_type":"markdown","source":"<h3> Dealing with Collinearity </h3>\n<a id=\"collinearity\"></a>\n<b> What we will do: </b>\n<ul>\n    <li><b>Remove the stress level feature</b> These two features are highly correlated, doesen't suprise me because stress levels is directly derived from the amount of children a patient has.  Since children was in our original dataframe, we should delete the stress level feature and I am still not sure the stress level feature provides additional insight to the dataset (Still working on my feature engineering skills).</li>\n    <li> <b> Adding features manually: </b> We will add features manually into our linear model and see how well it performs, also this will allow us to avoid collinearity </i>\n    <li><b>Repetitive Process: </b> I know it is a repetitive process nevertheless, my intent is to show how to avoid collinearity in our model step by step.</li>\n    </ul>\n\nKey: The main issue was the creation of the stress level variable which by removing that feature we can guarantee that collinearity is removed from our model."},{"metadata":{"trusted":true,"_uuid":"c46bb14a682cfbb04f1d8e8e601537dc9f75bcb8","_kg_hide-input":true},"cell_type":"code","source":"# Drop the columns we created.\n\ndf = df.drop(['weight_condition', 'charge_status', 'age_cat'], axis=1)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"1014a532a86719d999e98d9b5f82a123aa80fc28","_kg_hide-input":true},"cell_type":"code","source":"# Drop the columns that are highly correlated.\n\ncorr = df.corr()\ncorr","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"4bc32c3acafc740d7d360ebf6c6fbe84e90fe3a8","_kg_hide-input":true},"cell_type":"code","source":"# Drop stress_level\ndf.drop(['stress_level'], axis=1, inplace=True)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"c11cb2287fb4d26ba6c7584c03b983ff1dbd34b4","_kg_hide-input":true},"cell_type":"code","source":"fig = plt.figure(figsize=(12,8))\n\ncorr = df.corr()\nax = sns.heatmap(corr, linewidths=.5, cmap=\"RdBu\", annot=True, fmt=\"g\")\nplt.title(\"Correlation Plot \\n to Decide which features to drop\", fontsize=16)\nplt.show()","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"9c28703b8bd99560250f4a62e0ef99d0a58c1f8b"},"cell_type":"markdown","source":"<h2 align=\"center\">The Importance of Adding Smoking Status to our Model </h2>\n<a id=\"smoking_status\"></a>\n<h3>Model #1: Without the Smoking Feature: </h3>\n<ul>\n    <li><b>A low Adjusted R Squared:</b> We only have an adjusted R Square of 0.116, basically our model is not capturing any trend. </li>\n    <li><b>Components of our Model:</b> So only adding the Age and BMI features  to predict how much a Patient will be charged was not enough! </li>\n    <li><b>The importance of Exploratory Data Analysis:</b> Remember that in our exploratory analysis we didn't find any significant patterns with this features, so we should expect that model 1 will not perform as we would like to. But what if we add the smoking status? Remember, this feature was extremely useful in determining specific clusters in our sample population.  </li>\n</ul>\n\n<h3>Model #2: Adding the Smoking Feature: </h3>\n<ul> \n    <li><b>A higher Adjusted R Squared: </b> Just by adding the smoking feature we have an adjusted r squared of 0.747 as opposed to 0.116 without the smoking feature. </li>\n    <li><b>Adding the Smoking Feature: </b> Adding the smoking feature permits the model to have a better understanding of certain groups in our sample population. </li>\n    </ul>\n\n\n<b> Note: </b> Will give it some format throughout the weekend just wanted to show you two things: One how to avoid collinearity and two the importance of the smoker status of the patient towards determining how much they will be charged. All the other technicalities regarding the summary of the model I publish it in future kernel updates. "},{"metadata":{"trusted":true,"_uuid":"f06b375bae453fcc52cb557a72293fe70900dced","_kg_hide-input":true},"cell_type":"code","source":"model_without_smoking = ols(\"charges ~ bmi + age\", data=train).fit()\nprint(model_without_smoking.summary())","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"333adabfd85b0a2f995fc0d5a4a69cbbe2f880dd","_kg_hide-input":true},"cell_type":"code","source":"# Age in out X-axis since it has a higher correlation with charges.\nfig = plt.figure(figsize=(12,8))\nfig = sm.graphics.plot_regress_exog(model_without_smoking, \"age\", fig=fig)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"5f6df835ce93f347899851e5c58995fd16f950f8","_kg_hide-input":true},"cell_type":"code","source":"fig, ax = plt.subplots(figsize=(12, 8))\nfig = sm.graphics.plot_fit(model_without_smoking, \"age\", ax=ax)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"56cc9fc948729062c10abf4219cd96ea7548de39","_kg_hide-input":true},"cell_type":"code","source":"model_with_smoking = ols(\"charges ~ smoker + bmi + age\", data=train).fit()\nprint(model_with_smoking.summary())","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"e2237a6b7baa2d637683f1ae7e69e4b01ca18678","_kg_hide-input":true},"cell_type":"code","source":"# Age in out X-axis since it has a higher correlation with charges.\nfig = plt.figure(figsize=(12,8))\nfig = sm.graphics.plot_regress_exog(model_with_smoking, \"age\", fig=fig)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"72a43d9224a735356c08b1b6b7f7c74f846a7576","_kg_hide-input":true},"cell_type":"code","source":"fig, ax = plt.subplots(figsize=(12, 8))\nfig = sm.graphics.plot_fit(model_with_smoking, \"age\", ax=ax)","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"78aba4506e8853495dcc26f45a1fa441f5c757f9"},"cell_type":"markdown","source":"<h2> Dealing with Skewness </h2>\n<ul>\n    <li><b>Right-Skewness: </b> If the value is significantly greater than 0, the distribution should be more inclined to be more right skewed.</li>\n    <li><b>Symmetric: </b>  The closer the value of our skew is to 0, the more symetric is our distribution.</li>\n    <li><b>Left-Skewness:</b> If the value is significantly less than 0, the distribution should be more inclined to be left-skewed. </li>\n    </ul>\n    \n    \n <h3> How to find Skewness in our Distribution: </h3>   \n Coefficient of Skewness =$\\LARGE \\frac{\\overline{x} - z}{\\sigma}$<br><br>\n Where: <br>\n $\\overline{x}$ = Sample Mean <br>\n z = Mode <br>\n $\\sigma$ = Standard Deviation <br><br>\n \n <b> Standard Deviation: </b><br>\n $\\LARGE\\sigma{} = \\sqrt{\\frac{\\sum{x_i - \\overline{x}}}{n - 1}}$ <br>\n Where:<br>\n $\\large x_i$ = Observation <br>\n $\\large \\overline{x}$ = Sample Mean <br>\n $\\large n $ = Number of observations <br><br>\n \n <h3> Symmetric Distribution and Dealing with Outliers </h3>\nThere are many reasons why we would prefer to have a scaled label but one of the main reasons is to scaled down some observations that might be considered outliers. I will further explain this topic but didn't feel like leaving this section empty. Will provide more specific details in the future."},{"metadata":{"trusted":true,"_uuid":"2a8fbe348604512cb7e3cd79ae4447db9a8fb449","_kg_hide-input":true},"cell_type":"code","source":"# Let's see the skewness of charges\nnot_normalized = skew(df['charges'].values.tolist())\nnormalized = skew(np.log(df['charges'].values.tolist()))\n\n\ntrace0 = go.Bar(\n    x=['Not Normalized', 'Normalized'],\n    y=[not_normalized, normalized],\n    text=['Not Normalized Skewness', 'Normalized Skewness'],\n    marker=dict(\n        color='rgb(158,202,225)',\n        line=dict(\n            color='rgb(8,48,107)',\n            width=1.5,\n        )\n    ),\n    opacity=0.6\n)\n\ndata = [trace0]\nlayout = go.Layout(\n    title='Patient Charges Skewness \\n Normalized vs Not Normalized',\n    yaxis=dict(\n        title='Skewness Coeficient',\n        titlefont=dict(\n            size=16,\n            color='rgb(107, 107, 107)'\n        )\n))\n\nfig = go.Figure(data=data, layout=layout)\n\niplot(fig, filename='bar-direct-labels')","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"761710e6612563f7620d35d1678f94a81423d6b2","_kg_hide-input":true},"cell_type":"code","source":"plt.style.use('dark_background')\n\ncharges = df['charges'].values\n\nf, ((ax1, ax2), (ax3, ax4)) = plt.subplots(nrows=2, ncols=2, figsize=(20,18)) \n\n# Distribution of charges\nsns.distplot(train['charges'], ax=ax1)\nax1.set_title(\"Distribution of Charges \\n (Right-Skewed)\", fontsize=16)\n\n\n# Scaling charges using a natural logarithm.\nsns.distplot(np.log(train['charges']), ax=ax2)\nax2.set_title(\"Scaled Distribution of Charges \\n (Symmetric Distribution)\", fontsize=16)\n\n# Quantile plot \nstats.probplot(train[\"charges\"], plot=ax3)\nax3.set_title(\"Probability Plot \\n (Unscaled)\", fontsize=16)\n\nstats.probplot(np.log(train[\"charges\"]), plot=ax4)\nax4.set_title(\"Probability Plot \\n (Scaled using Natural Logarithm)\", fontsize=16)\n\nplt.show()","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"a936957e3149d4a1469dc9ea54906e29742b5188"},"cell_type":"markdown","source":"### Using Our Model with a Symmetric Patient Charge's Distribution:\n---> To be Updated"},{"metadata":{"_kg_hide-input":true,"trusted":true,"_uuid":"eb3f622468b263f1f779a469b97cfd4b304e547a"},"cell_type":"code","source":"train['log_charges'] = np.log(train[\"charges\"])\n\nmodel_with_logcharges = ols(\"log_charges ~ smoker + bmi + age\", data=train).fit()\nprint(model_with_smoking.summary())","execution_count":null,"outputs":[]},{"metadata":{"_kg_hide-input":true,"trusted":true,"_uuid":"445e56b76ebcef023d2a2db2d831e1c6f702c166"},"cell_type":"code","source":"# Using age to predict charges\nplt.style.use(\"dark_background\")\n\nfig = plt.figure(figsize=(12,8))\nfig = sm.graphics.plot_regress_exog(model_with_logcharges, \"age\", fig=fig)","execution_count":null,"outputs":[]},{"metadata":{"_kg_hide-input":true,"trusted":true,"_uuid":"c1974e66b4c43311529b629f43950daa501a2bf1"},"cell_type":"code","source":"fig, ax = plt.subplots(figsize=(12, 8))\nfig = sm.graphics.plot_fit(model_with_logcharges, \"age\", ax=ax)","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"1747fb7d1c6698bb105dc465d83f03ee7a952c9d"},"cell_type":"markdown","source":"### Feature Engineering with Itertools (Up Next):"}],"metadata":{"kernelspec":{"display_name":"Python 3","language":"python","name":"python3"},"language_info":{"name":"python","version":"3.6.6","mimetype":"text/x-python","codemirror_mode":{"name":"ipython","version":3},"pygments_lexer":"ipython3","nbconvert_exporter":"python","file_extension":".py"}},"nbformat":4,"nbformat_minor":1}