Multi Linear Regression (1).ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import pandas as pd\n",
    "import matplotlib.pyplot as plt\n",
    "import seaborn as sns\n",
    "from statsmodels.graphics.regressionplots import influence_plot\n",
    "import statsmodels.formula.api as smf\n",
    "import numpy as np"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Read the data\n",
    "cars = pd.read_csv(\"Cars.csv\")\n",
    "cars.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "cars.info()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#check for missing values\n",
    "cars.isna().sum()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Correlation Matrix"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "cars.corr()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Scatterplot between variables along with histograms"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Format the plot background and scatter plots for all the variables\n",
    "sns.set_style(style='darkgrid')\n",
    "sns.pairplot(cars)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Preparing a model"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Build model\n",
    "import statsmodels.formula.api as smf \n",
    "model = smf.ols('MPG~WT+VOL+SP+HP',data=cars).fit()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Coefficients\n",
    "model.params"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#t and p-Values\n",
    "print(model.tvalues, '\\n', model.pvalues)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#R squared values\n",
    "(model.rsquared,model.rsquared_adj)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Simple Linear Regression Models"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "ml_v=smf.ols('MPG~VOL',data = cars).fit()  \n",
    "#t and p-Values\n",
    "print(ml_v.tvalues, '\\n', ml_v.pvalues)  "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "ml_w=smf.ols('MPG~WT',data = cars).fit()  \n",
    "print(ml_w.tvalues, '\\n', ml_w.pvalues)  "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "ml_wv=smf.ols('MPG~WT+VOL',data = cars).fit()  \n",
    "print(ml_wv.tvalues, '\\n', ml_wv.pvalues)  "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Calculating VIF"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "rsq_hp = smf.ols('HP~WT+VOL+SP',data=cars).fit().rsquared  \n",
    "vif_hp = 1/(1-rsq_hp) # 16.33\n",
    "\n",
    "rsq_wt = smf.ols('WT~HP+VOL+SP',data=cars).fit().rsquared  \n",
    "vif_wt = 1/(1-rsq_wt) # 564.98\n",
    "\n",
    "rsq_vol = smf.ols('VOL~WT+SP+HP',data=cars).fit().rsquared  \n",
    "vif_vol = 1/(1-rsq_vol) #  564.84\n",
    "\n",
    "rsq_sp = smf.ols('SP~WT+VOL+HP',data=cars).fit().rsquared  \n",
    "vif_sp = 1/(1-rsq_sp) #  16.35\n",
    "\n",
    "# Storing vif values in a data frame\n",
    "d1 = {'Variables':['Hp','WT','VOL','SP'],'VIF':[vif_hp,vif_wt,vif_vol,vif_sp]}\n",
    "Vif_frame = pd.DataFrame(d1)  \n",
    "Vif_frame"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Residual Analysis"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Test for Normality of Residuals (Q-Q Plot)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import statsmodels.api as sm\n",
    "qqplot=sm.qqplot(model.resid,line='q') # line = 45 to draw the diagnoal line\n",
    "plt.title(\"Normal Q-Q plot of residuals\")\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "list(np.where(model.resid>10))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Residual Plot for Homoscedasticity"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def get_standardized_values( vals ):\n",
    "    return (vals - vals.mean())/vals.std()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plt.scatter(get_standardized_values(model.fittedvalues),\n",
    "            get_standardized_values(model.resid))\n",
    "\n",
    "plt.title('Residual Plot')\n",
    "plt.xlabel('Standardized Fitted values')\n",
    "plt.ylabel('Standardized residual values')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Residual Vs Regressors"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "fig = plt.figure(figsize=(15,8))\n",
    "fig = sm.graphics.plot_regress_exog(model, \"VOL\", fig=fig)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "fig = plt.figure(figsize=(15,8))\n",
    "fig = sm.graphics.plot_regress_exog(model, \"SP\", fig=fig)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "fig = plt.figure(figsize=(15,8))\n",
    "fig = sm.graphics.plot_regress_exog(model, \"HP\", fig=fig)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "fig = plt.figure(figsize=(15,8))\n",
    "fig = sm.graphics.plot_regress_exog(model, \"WT\", fig=fig)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Model Deletion Diagnostics"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Detecting Influencers/Outliers"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Cook’s Distance"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "model_influence = model.get_influence()\n",
    "(c,_ ) = model_influence.cooks_distance"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "c"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Plot the influencers values using stem plot\n",
    "fig = plt.subplots(figsize=(20, 7))\n",
    "plt.stem(np.arange(len(cars)), np.round(c, 3))\n",
    "plt.xlabel('Row index')\n",
    "plt.ylabel('Cooks Distance')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#index and value of influencer where c is more than .5\n",
    "(np.argmax(c),np.max(c))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## High Influence points"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from statsmodels.graphics.regressionplots import influence_plot\n",
    "influence_plot(model)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "k = cars.shape[1]\n",
    "n = cars.shape[0]\n",
    "leverage_cutoff = 3*((k + 1)/n)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "leverage_cutoff"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### From the above plot, it is evident that data point 70 and 76 are the influencers"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "cars[cars.index.isin([70, 76])]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#See the differences in HP and other variable values\n",
    "cars.head()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Improving the model"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Load the data\n",
    "cars_new = pd.read_csv(\"Cars.csv\",index_col=0)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Discard the data points which are influencers and reasign the row number (reset_index())\n",
    "car1=cars_new.drop(cars_new.index[[70,76]],axis=0).reset_index()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "car1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Build Model"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Exclude variable \"WT\" and generate R-Squared and AIC values\n",
    "final_ml_V= smf.ols('MPG~VOL+SP+HP',data = car1).fit()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "(final_ml_V.rsquared,final_ml_V.aic)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Exclude variable \"VOL\" and generate R-Squared and AIC values\n",
    "final_ml_W= smf.ols('MPG~WT+SP+HP',data = car1).fit()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "(final_ml_W.rsquared,final_ml_W.aic)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### Comparing above R-Square and AIC values, model 'final_ml_V' has high R- square and low AIC value hence include variable 'VOL' so that multi collinearity problem would be resolved."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Cook’s Distance"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "model_influence_V = final_ml_V.get_influence()\n",
    "(c_V, _) = model_influence_V.cooks_distance"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "fig= plt.subplots(figsize=(20,7))\n",
    "plt.stem(np.arange(len(car1)),np.round(c_V,3));\n",
    "plt.xlabel('Row index')\n",
    "plt.ylabel('Cooks Distance');"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#index of the data points where c is more than .5\n",
    "(np.argmax(c_V),np.max(c_V))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Drop 76 and 77 observations\n",
    "car2=car1.drop(car1.index[[76,77]],axis=0)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "car2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Reset the index and re arrange the row values\n",
    "car3=car2.reset_index()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "car4=car3.drop(['index'],axis=1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "car4"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Build the model on the new data\n",
    "final_ml_V= smf.ols('MPG~VOL+SP+HP',data = car4).fit()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Again check for influencers\n",
    "model_influence_V = final_ml_V.get_influence()\n",
    "(c_V, _) = model_influence_V.cooks_distance"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "fig= plt.subplots(figsize=(20,7))\n",
    "plt.stem(np.arange(len(car4)),np.round(c_V,3));\n",
    "plt.xlabel('Row index')\n",
    "plt.ylabel('Cooks Distance');"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#index of the data points where c is more than .5\n",
    "(np.argmax(c_V),np.max(c_V))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Since the value is <1 , we can stop the diagnostic process and finalize the model"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Check the accuracy of the mode\n",
    "final_ml_V= smf.ols('MPG~VOL+SP+HP',data = car4).fit()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "(final_ml_V.rsquared,final_ml_V.aic)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Predicting for new data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#New data for prediction\n",
    "new_data=pd.DataFrame({'HP':40,\"VOL\":95,\"SP\":102},index=[1])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "new_data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "final_ml_V.predict(new_data)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "pred_y = final_ml_V.predict(car4)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "pred_y"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "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.7.4"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}