{"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"pygments_lexer":"ipython3","nbconvert_exporter":"python","version":"3.6.4","file_extension":".py","codemirror_mode":{"name":"ipython","version":3},"name":"python","mimetype":"text/x-python"}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"code","source":"!pip install probscale","metadata":{"execution":{"iopub.status.busy":"2022-08-11T17:38:49.982725Z","iopub.execute_input":"2022-08-11T17:38:49.983278Z","iopub.status.idle":"2022-08-11T17:39:04.785423Z","shell.execute_reply.started":"2022-08-11T17:38:49.983178Z","shell.execute_reply":"2022-08-11T17:39:04.784151Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Import the necessary libraries\nimport csv\nimport re\nimport numpy as np\nimport seaborn as sns\nimport matplotlib.pyplot as plt\nimport pandas as pd\nimport probscale\n\nplt.style.use('ggplot')\n\n%config InlineBackend.figure_format = 'retina'\n%matplotlib inline","metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","execution":{"iopub.status.busy":"2022-08-11T17:39:04.788659Z","iopub.execute_input":"2022-08-11T17:39:04.789740Z","iopub.status.idle":"2022-08-11T17:39:05.499971Z","shell.execute_reply.started":"2022-08-11T17:39:04.789679Z","shell.execute_reply":"2022-08-11T17:39:05.499076Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import statsmodels.api as sm\nfrom sklearn import metrics, linear_model\nfrom sklearn.model_selection import KFold, train_test_split, cross_val_score\nfrom sklearn.preprocessing import StandardScaler\nfrom sklearn.linear_model import LinearRegression, Ridge, Lasso, RidgeCV, LassoCV\nfrom statsmodels.graphics.gofplots import qqplot\nfrom scipy import stats\nfrom scipy.stats import shapiro","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:57:42.336772Z","iopub.execute_input":"2022-08-11T18:57:42.337404Z","iopub.status.idle":"2022-08-11T18:57:43.079786Z","shell.execute_reply.started":"2022-08-11T18:57:42.337356Z","shell.execute_reply":"2022-08-11T18:57:43.078868Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Read Data into DataFrame","metadata":{}},{"cell_type":"code","source":"df = pd.read_csv('../input/house-prices-advanced-regression-techniques/train.csv')\ndf","metadata":{"execution":{"iopub.status.busy":"2022-08-11T17:39:05.692880Z","iopub.execute_input":"2022-08-11T17:39:05.694088Z","iopub.status.idle":"2022-08-11T17:39:05.774017Z","shell.execute_reply.started":"2022-08-11T17:39:05.694036Z","shell.execute_reply":"2022-08-11T17:39:05.772832Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# DataFrame contains lots of blank entries that needs to be turned into proper NaN\ndf = df.replace(['NA', ''], np.NaN)\ndf","metadata":{"execution":{"iopub.status.busy":"2022-08-11T17:39:05.775828Z","iopub.execute_input":"2022-08-11T17:39:05.776598Z","iopub.status.idle":"2022-08-11T17:39:05.820499Z","shell.execute_reply.started":"2022-08-11T17:39:05.776551Z","shell.execute_reply":"2022-08-11T17:39:05.819373Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"https://towardsdatascience.com/wrangling-through-dataland-modeling-house-prices-in-ames-iowa-75b9b4086c96","metadata":{}},{"cell_type":"code","source":"# Identify the numeric data columns for the data description file\ncols = ['MSSubClass', 'LotFrontage', 'LotArea', 'OverallQual', 'OverallCond', 'YearBuilt', 'YearRemodAdd', \n        'MasVnrArea', 'BsmtFinSF1', 'BsmtFinSF2', 'BsmtUnfSF', 'TotalBsmtSF', '1stFlrSF', '2ndFlrSF', \n        'LowQualFinSF', 'GrLivArea', 'BsmtFullBath', 'BsmtHalfBath', 'FullBath', 'HalfBath', 'BedroomAbvGr', \n        'KitchenAbvGr', 'TotRmsAbvGrd', 'Fireplaces', 'GarageYrBlt', 'GarageCars', 'GarageArea', 'WoodDeckSF', \n        'OpenPorchSF', 'EnclosedPorch', '3SsnPorch', 'ScreenPorch', 'PoolArea', 'MiscVal', 'MoSold', 'YrSold', \n        'SalePrice']","metadata":{"execution":{"iopub.status.busy":"2022-08-11T17:39:05.822455Z","iopub.execute_input":"2022-08-11T17:39:05.822784Z","iopub.status.idle":"2022-08-11T17:39:05.828761Z","shell.execute_reply.started":"2022-08-11T17:39:05.822755Z","shell.execute_reply":"2022-08-11T17:39:05.827912Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Preview the list of the columns\ndf.columns","metadata":{"execution":{"iopub.status.busy":"2022-08-11T17:39:05.830235Z","iopub.execute_input":"2022-08-11T17:39:05.830642Z","iopub.status.idle":"2022-08-11T17:39:05.843644Z","shell.execute_reply.started":"2022-08-11T17:39:05.830608Z","shell.execute_reply":"2022-08-11T17:39:05.842841Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Cast the numeric columns within the dataframe to a float datatype\ndf[cols] = df[cols].astype('float')","metadata":{"execution":{"iopub.status.busy":"2022-08-11T17:39:05.845111Z","iopub.execute_input":"2022-08-11T17:39:05.845643Z","iopub.status.idle":"2022-08-11T17:39:05.871039Z","shell.execute_reply.started":"2022-08-11T17:39:05.845612Z","shell.execute_reply":"2022-08-11T17:39:05.870093Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Preview the info of the columns with their datatypes and non-null content\ndf.info()","metadata":{"execution":{"iopub.status.busy":"2022-08-11T17:39:05.872292Z","iopub.execute_input":"2022-08-11T17:39:05.872992Z","iopub.status.idle":"2022-08-11T17:39:05.905113Z","shell.execute_reply.started":"2022-08-11T17:39:05.872928Z","shell.execute_reply":"2022-08-11T17:39:05.903792Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Print statistics of the dataframe\ndf.describe()","metadata":{"execution":{"iopub.status.busy":"2022-08-11T17:39:05.909614Z","iopub.execute_input":"2022-08-11T17:39:05.910331Z","iopub.status.idle":"2022-08-11T17:39:06.022487Z","shell.execute_reply.started":"2022-08-11T17:39:05.910282Z","shell.execute_reply":"2022-08-11T17:39:06.021129Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Data Engineering","metadata":{}},{"cell_type":"code","source":"# Dropping those variables with less than 300 non-null values\ndf.drop(['Alley', 'FireplaceQu', 'PoolQC', 'Fence', 'MiscFeature'], axis=1, inplace=True)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T17:39:06.025024Z","iopub.execute_input":"2022-08-11T17:39:06.025948Z","iopub.status.idle":"2022-08-11T17:39:06.034959Z","shell.execute_reply.started":"2022-08-11T17:39:06.025894Z","shell.execute_reply":"2022-08-11T17:39:06.033506Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Preview the number of rows and columns within the dataframe\ndf.shape","metadata":{"execution":{"iopub.status.busy":"2022-08-11T17:39:06.037012Z","iopub.execute_input":"2022-08-11T17:39:06.037796Z","iopub.status.idle":"2022-08-11T17:39:06.048669Z","shell.execute_reply.started":"2022-08-11T17:39:06.037718Z","shell.execute_reply":"2022-08-11T17:39:06.047164Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df.MSZoning.unique() # C denotes commercial properties","metadata":{"execution":{"iopub.status.busy":"2022-08-11T17:39:06.050757Z","iopub.execute_input":"2022-08-11T17:39:06.051133Z","iopub.status.idle":"2022-08-11T17:39:06.060195Z","shell.execute_reply.started":"2022-08-11T17:39:06.051101Z","shell.execute_reply":"2022-08-11T17:39:06.059084Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Delete all data with MSZoning = commercial, agriculture and industrial as these are not residential units\ndf = df[(df.MSZoning != 'C (all)') & (df.MSZoning != 'I (all)') & (df.MSZoning != 'A (agr)')]","metadata":{"execution":{"iopub.status.busy":"2022-08-11T17:39:06.061503Z","iopub.execute_input":"2022-08-11T17:39:06.062694Z","iopub.status.idle":"2022-08-11T17:39:06.073545Z","shell.execute_reply.started":"2022-08-11T17:39:06.062658Z","shell.execute_reply":"2022-08-11T17:39:06.072677Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Notice a large number of null values in 'LotFrontage'. Examine the 'Lots' grouping.\ndf_lots = df[['LotFrontage', 'LotArea', 'LotConfig', 'LotShape']]\ngrouped_lots = df_lots.groupby(['LotShape'])\ngrouped_lots.mean()","metadata":{"execution":{"iopub.status.busy":"2022-08-11T17:39:06.075613Z","iopub.execute_input":"2022-08-11T17:39:06.076006Z","iopub.status.idle":"2022-08-11T17:39:06.099395Z","shell.execute_reply.started":"2022-08-11T17:39:06.075973Z","shell.execute_reply":"2022-08-11T17:39:06.098033Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"grouped_lots2 = df_lots.groupby(['LotConfig'])\ngrouped_lots2.count()","metadata":{"execution":{"iopub.status.busy":"2022-08-11T17:39:06.101372Z","iopub.execute_input":"2022-08-11T17:39:06.101736Z","iopub.status.idle":"2022-08-11T17:39:06.116245Z","shell.execute_reply.started":"2022-08-11T17:39:06.101704Z","shell.execute_reply":"2022-08-11T17:39:06.115356Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Generate 'Lots' group where there are null 'LotFrontage' values\ndf_LotFrontage_NA = df_lots.loc[(df['LotFrontage'].isnull())]\ndf_LotFrontage_NA.columns","metadata":{"execution":{"iopub.status.busy":"2022-08-11T17:39:06.117480Z","iopub.execute_input":"2022-08-11T17:39:06.118592Z","iopub.status.idle":"2022-08-11T17:39:06.126841Z","shell.execute_reply.started":"2022-08-11T17:39:06.118555Z","shell.execute_reply":"2022-08-11T17:39:06.125667Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Count the number of rows containing null LotFrontage values\ndf_LotFrontage_NA.LotFrontage.isnull().sum()","metadata":{"execution":{"iopub.status.busy":"2022-08-11T17:39:06.128668Z","iopub.execute_input":"2022-08-11T17:39:06.129583Z","iopub.status.idle":"2022-08-11T17:39:06.140522Z","shell.execute_reply.started":"2022-08-11T17:39:06.129535Z","shell.execute_reply":"2022-08-11T17:39:06.139352Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# A reasonable assumption is that LotFrontage is linked to LotConfig and LotShape, and\n# the other 'Lot' variables have all observations\n# So all NA in 'LotFrontage' will be replaced with its mean based on the corresponding 'LotShape', which is indexed at 3\ndf_LotFrontage_NA['LotFrontage'] =  df_LotFrontage_NA.apply(lambda x: 74.7688 if (x[3] == 'IR1') \n                                                           else x[0], axis=1)\ndf_LotFrontage_NA['LotFrontage'] = df_LotFrontage_NA.apply(lambda x: 67.4375 if (x[3] == 'IR2') \n                                                           else x[0], axis=1)\ndf_LotFrontage_NA['LotFrontage'] = df_LotFrontage_NA.apply(lambda x: 117.6364 if (x[3] == 'IR3') \n                                                           else x[0], axis=1)\ndf_LotFrontage_NA['LotFrontage'] = df_LotFrontage_NA.apply(lambda x: 66.8214 if (x[3] == 'Reg') \n                                                           else x[0], axis=1)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T17:39:06.141826Z","iopub.execute_input":"2022-08-11T17:39:06.142758Z","iopub.status.idle":"2022-08-11T17:39:06.168696Z","shell.execute_reply.started":"2022-08-11T17:39:06.142721Z","shell.execute_reply":"2022-08-11T17:39:06.167859Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Count the LotFrontage with NA values\nLotFront_fills = df_LotFrontage_NA.LotFrontage\nLotFront_fills.count()","metadata":{"execution":{"iopub.status.busy":"2022-08-11T17:39:06.169950Z","iopub.execute_input":"2022-08-11T17:39:06.170551Z","iopub.status.idle":"2022-08-11T17:39:06.177288Z","shell.execute_reply.started":"2022-08-11T17:39:06.170516Z","shell.execute_reply":"2022-08-11T17:39:06.176362Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"LotFront_fills","metadata":{"execution":{"iopub.status.busy":"2022-08-11T17:44:29.346277Z","iopub.execute_input":"2022-08-11T17:44:29.346735Z","iopub.status.idle":"2022-08-11T17:44:29.356902Z","shell.execute_reply.started":"2022-08-11T17:44:29.346701Z","shell.execute_reply":"2022-08-11T17:44:29.355624Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Sum the NaN values of the LotFrontage in the original dataset to verify\ndf.LotFrontage.isnull().sum()","metadata":{"execution":{"iopub.status.busy":"2022-08-11T17:39:06.178417Z","iopub.execute_input":"2022-08-11T17:39:06.179383Z","iopub.status.idle":"2022-08-11T17:39:06.190335Z","shell.execute_reply.started":"2022-08-11T17:39:06.179336Z","shell.execute_reply":"2022-08-11T17:39:06.189001Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Fill the 'LotFrontage' null values with the LotFront_fills series in the given order of the data\ndf.loc[df.LotFrontage.isnull(), 'LotFrontage'] = LotFront_fills\ndf.LotFrontage.isnull().sum()","metadata":{"execution":{"iopub.status.busy":"2022-08-11T17:45:57.613648Z","iopub.execute_input":"2022-08-11T17:45:57.614095Z","iopub.status.idle":"2022-08-11T17:45:57.625802Z","shell.execute_reply.started":"2022-08-11T17:45:57.614058Z","shell.execute_reply":"2022-08-11T17:45:57.624673Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Adding structure age variable\ndf['Age'] = df.apply(lambda x: x['YrSold']-x['YearBuilt'] if (x['YearBuilt']<x['YearRemodAdd']) \n                                                           else (x['YrSold']-x['YearRemodAdd']), axis=1)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T17:47:22.412081Z","iopub.execute_input":"2022-08-11T17:47:22.412594Z","iopub.status.idle":"2022-08-11T17:47:22.460807Z","shell.execute_reply.started":"2022-08-11T17:47:22.412554Z","shell.execute_reply":"2022-08-11T17:47:22.459657Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# The right-skew to the SalePrice is obvious in the histogram, showing that SalePrice has a long right tail\nfig, ax = plt.subplots(figsize=(10,5))\n\nsns.distplot(df.SalePrice, bins=30, kde=True, ax=ax, color='forestgreen')\nplt.title('Housing Sale Price Histogram', fontsize=15)\nplt.xlabel('Sale Price, $', fontsize=12)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T17:48:04.345689Z","iopub.execute_input":"2022-08-11T17:48:04.346975Z","iopub.status.idle":"2022-08-11T17:48:04.805062Z","shell.execute_reply.started":"2022-08-11T17:48:04.346923Z","shell.execute_reply":"2022-08-11T17:48:04.803823Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Creating variables on property area/size","metadata":{}},{"cell_type":"code","source":"# One of the chief factors affecting house prices in the size of the property.\n# Here we examine the basement square footage variables.\ndf[['BsmtFinSF1', 'BsmtFinSF2', 'TotalBsmtSF', 'BsmtUnfSF']].head()","metadata":{"execution":{"iopub.status.busy":"2022-08-11T17:52:24.209450Z","iopub.execute_input":"2022-08-11T17:52:24.209978Z","iopub.status.idle":"2022-08-11T17:52:24.228124Z","shell.execute_reply.started":"2022-08-11T17:52:24.209938Z","shell.execute_reply":"2022-08-11T17:52:24.226945Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Here we examine the ground square footage variables. The first column is obviously the sum of the other two.\ndf[['GrLivArea', '1stFlrSF', '2ndFlrSF']].head()","metadata":{"execution":{"iopub.status.busy":"2022-08-11T17:53:41.983296Z","iopub.execute_input":"2022-08-11T17:53:41.983789Z","iopub.status.idle":"2022-08-11T17:53:41.998950Z","shell.execute_reply.started":"2022-08-11T17:53:41.983748Z","shell.execute_reply":"2022-08-11T17:53:41.997588Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df.GrLivArea.describe()","metadata":{"execution":{"iopub.status.busy":"2022-08-11T17:54:09.603451Z","iopub.execute_input":"2022-08-11T17:54:09.603901Z","iopub.status.idle":"2022-08-11T17:54:09.617133Z","shell.execute_reply.started":"2022-08-11T17:54:09.603845Z","shell.execute_reply":"2022-08-11T17:54:09.615990Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Create a total \"finished\" basement square footage variable\ndf['BaseLivArea'] = df.TotalBsmtSF - df.BsmtUnfSF","metadata":{"execution":{"iopub.status.busy":"2022-08-11T17:56:56.783057Z","iopub.execute_input":"2022-08-11T17:56:56.783533Z","iopub.status.idle":"2022-08-11T17:56:56.791509Z","shell.execute_reply.started":"2022-08-11T17:56:56.783494Z","shell.execute_reply":"2022-08-11T17:56:56.790677Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Clear positive relationship between 'TSF' and 'SalePrice' variables, which is what one would expect.\n# The larger the house, the higher its price, ceteris paribus\nfig, ax = plt.subplots(figsize=(10,5))\n\nsns.regplot(x=\"GrLivArea\", y=\"SalePrice\", data=df, ax=ax)\nax.set_title('Sale Price vs Above Grade Square Footage')","metadata":{"execution":{"iopub.status.busy":"2022-08-11T17:57:59.423678Z","iopub.execute_input":"2022-08-11T17:57:59.424129Z","iopub.status.idle":"2022-08-11T17:57:59.968967Z","shell.execute_reply.started":"2022-08-11T17:57:59.424095Z","shell.execute_reply":"2022-08-11T17:57:59.967716Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df['BaseLivArea']","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:00:21.675965Z","iopub.execute_input":"2022-08-11T18:00:21.677198Z","iopub.status.idle":"2022-08-11T18:00:21.686982Z","shell.execute_reply.started":"2022-08-11T18:00:21.677151Z","shell.execute_reply":"2022-08-11T18:00:21.685990Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# A simple linear regression of TSF on the SalePrice gives an R-squared of 50% by itself!\nols = LinearRegression()\n\nX_GLA = df[['GrLivArea']]\ny_SalePrice = df.SalePrice\n\nols.fit(X_GLA, y_SalePrice)\nprint(\"R-squared of Total Square Footage on Sale Price:\", ols.score(X_GLA, y_SalePrice))","metadata":{"execution":{"iopub.status.busy":"2022-08-11T17:59:18.193726Z","iopub.execute_input":"2022-08-11T17:59:18.194199Z","iopub.status.idle":"2022-08-11T17:59:18.216652Z","shell.execute_reply.started":"2022-08-11T17:59:18.194162Z","shell.execute_reply":"2022-08-11T17:59:18.215578Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# There also appears to be a positive relationship between 'GrLivArea' and 'BaseLivArea', plus\n# there seems to be greater variation in 'BaseLivArea' as 'GrLivArea' increases\nfig, ax = plt.subplots(figsize=(10,5))\n\nsns.regplot(x=\"GrLivArea\", y=\"BaseLivArea\", data=df, ax=ax) \nax.set_title('Above Grade Square Footage vs Basement Square Footage')","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:00:56.911672Z","iopub.execute_input":"2022-08-11T18:00:56.912097Z","iopub.status.idle":"2022-08-11T18:00:57.417099Z","shell.execute_reply.started":"2022-08-11T18:00:56.912063Z","shell.execute_reply":"2022-08-11T18:00:57.416016Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df[['SalePrice', 'GrLivArea', 'BaseLivArea']].describe()","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:14:31.355580Z","iopub.execute_input":"2022-08-11T18:14:31.356067Z","iopub.status.idle":"2022-08-11T18:14:31.382290Z","shell.execute_reply.started":"2022-08-11T18:14:31.356028Z","shell.execute_reply":"2022-08-11T18:14:31.381469Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Searching for outliers","metadata":{}},{"cell_type":"code","source":"# There are 22 instances of 'SalePrice' being 3 std above the mean\ndf.SalePrice[(df.SalePrice > np.mean(df.SalePrice) + 3*np.std(df.SalePrice))].count()","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:15:21.686642Z","iopub.execute_input":"2022-08-11T18:15:21.687192Z","iopub.status.idle":"2022-08-11T18:15:21.700446Z","shell.execute_reply.started":"2022-08-11T18:15:21.687150Z","shell.execute_reply":"2022-08-11T18:15:21.699259Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# There are zero cases of the opposite\ndf.SalePrice[(df.SalePrice < np.mean(df.SalePrice) - 3*np.std(df.SalePrice))].count()","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:16:06.374678Z","iopub.execute_input":"2022-08-11T18:16:06.375553Z","iopub.status.idle":"2022-08-11T18:16:06.386392Z","shell.execute_reply.started":"2022-08-11T18:16:06.375507Z","shell.execute_reply":"2022-08-11T18:16:06.384926Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# There are 16 instances of 'GrLivArea' being 3 std above the mean\ndf.GrLivArea[(df.GrLivArea > np.mean(df.GrLivArea) + 3*np.std(df.GrLivArea))].count()","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:16:37.034705Z","iopub.execute_input":"2022-08-11T18:16:37.035158Z","iopub.status.idle":"2022-08-11T18:16:37.044993Z","shell.execute_reply.started":"2022-08-11T18:16:37.035121Z","shell.execute_reply":"2022-08-11T18:16:37.043699Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# But zero of the opposite\ndf.GrLivArea[(df.GrLivArea < np.mean(df.GrLivArea) - 3*np.std(df.GrLivArea))].count()","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:17:20.614628Z","iopub.execute_input":"2022-08-11T18:17:20.615611Z","iopub.status.idle":"2022-08-11T18:17:20.625007Z","shell.execute_reply.started":"2022-08-11T18:17:20.615566Z","shell.execute_reply":"2022-08-11T18:17:20.623908Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# There are 6 instances of 'BaseLivArea' being 3 std above the mean\ndf.BaseLivArea[(df.BaseLivArea > np.mean(df.BaseLivArea) + 3*np.std(df.BaseLivArea))].count()","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:17:44.607098Z","iopub.execute_input":"2022-08-11T18:17:44.607487Z","iopub.status.idle":"2022-08-11T18:17:44.615596Z","shell.execute_reply.started":"2022-08-11T18:17:44.607456Z","shell.execute_reply":"2022-08-11T18:17:44.614463Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# But zero of the opposite\ndf.BaseLivArea[(df.BaseLivArea < np.mean(df.BaseLivArea) - 3*np.std(df.BaseLivArea))].count()","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:18:03.976517Z","iopub.execute_input":"2022-08-11T18:18:03.976952Z","iopub.status.idle":"2022-08-11T18:18:03.986811Z","shell.execute_reply.started":"2022-08-11T18:18:03.976916Z","shell.execute_reply":"2022-08-11T18:18:03.985417Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df.shape","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:18:35.774739Z","iopub.execute_input":"2022-08-11T18:18:35.775432Z","iopub.status.idle":"2022-08-11T18:18:35.782121Z","shell.execute_reply.started":"2022-08-11T18:18:35.775390Z","shell.execute_reply":"2022-08-11T18:18:35.780928Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_exout = df[(df.SalePrice < np.mean(df.SalePrice) + 3*np.std(df.SalePrice))]","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:18:59.356637Z","iopub.execute_input":"2022-08-11T18:18:59.357367Z","iopub.status.idle":"2022-08-11T18:18:59.369767Z","shell.execute_reply.started":"2022-08-11T18:18:59.357300Z","shell.execute_reply":"2022-08-11T18:18:59.368144Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_exout = df_exout[(df_exout.GrLivArea < np.mean(df.GrLivArea) + 3*np.std(df.GrLivArea))]","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:19:24.203996Z","iopub.execute_input":"2022-08-11T18:19:24.204470Z","iopub.status.idle":"2022-08-11T18:19:24.214550Z","shell.execute_reply.started":"2022-08-11T18:19:24.204434Z","shell.execute_reply":"2022-08-11T18:19:24.213650Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_exout = df_exout[(df_exout.BaseLivArea < np.mean(df.BaseLivArea) + 3*np.std(df.BaseLivArea))]","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:19:36.739246Z","iopub.execute_input":"2022-08-11T18:19:36.739955Z","iopub.status.idle":"2022-08-11T18:19:36.749548Z","shell.execute_reply.started":"2022-08-11T18:19:36.739906Z","shell.execute_reply":"2022-08-11T18:19:36.748584Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# All in all, 32 outliers are deleted\ndf_exout.shape ","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:19:57.841901Z","iopub.execute_input":"2022-08-11T18:19:57.842296Z","iopub.status.idle":"2022-08-11T18:19:57.849184Z","shell.execute_reply.started":"2022-08-11T18:19:57.842265Z","shell.execute_reply":"2022-08-11T18:19:57.848250Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_exout['PriceSF'] = df_exout.SalePrice / df_exout.GrLivArea","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:21:01.701571Z","iopub.execute_input":"2022-08-11T18:21:01.703036Z","iopub.status.idle":"2022-08-11T18:21:01.711164Z","shell.execute_reply.started":"2022-08-11T18:21:01.702977Z","shell.execute_reply":"2022-08-11T18:21:01.709953Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Searching for non-commercial transactions","metadata":{}},{"cell_type":"code","source":"print(df_exout.SalePrice.mean())\nprint(df_exout.PriceSF.mean())","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:21:30.895309Z","iopub.execute_input":"2022-08-11T18:21:30.895758Z","iopub.status.idle":"2022-08-11T18:21:30.904533Z","shell.execute_reply.started":"2022-08-11T18:21:30.895723Z","shell.execute_reply":"2022-08-11T18:21:30.903580Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"np.sum(df_exout.SaleCondition == 'Abnorml')","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:23:57.634708Z","iopub.execute_input":"2022-08-11T18:23:57.635153Z","iopub.status.idle":"2022-08-11T18:23:57.644922Z","shell.execute_reply.started":"2022-08-11T18:23:57.635115Z","shell.execute_reply":"2022-08-11T18:23:57.643724Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# 93 of these \"abnormal\" sale transactions are listed as \"abnormal\" or to family members,\n# and their mean SalePrice is well below the overall mean SalePrice\nprint(df_exout.SalePrice[(df_exout.SaleCondition == 'Abnorml') | (df_exout.SaleCondition == 'Family')].count())\ndf_exout.SalePrice[(df_exout.SaleCondition == 'Abnorml') | (df_exout.SaleCondition == 'Family')].mean()","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:26:28.314725Z","iopub.execute_input":"2022-08-11T18:26:28.316079Z","iopub.status.idle":"2022-08-11T18:26:28.330556Z","shell.execute_reply.started":"2022-08-11T18:26:28.316025Z","shell.execute_reply":"2022-08-11T18:26:28.329101Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Even the mean 'PriceSF' of these \"abnormal\" sales is well below the overall mean\ndf_exout.PriceSF[(df_exout.SaleCondition == 'Abnorml') | (df_exout.SaleCondition == 'Family')].mean()","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:27:08.450726Z","iopub.execute_input":"2022-08-11T18:27:08.451177Z","iopub.status.idle":"2022-08-11T18:27:08.462540Z","shell.execute_reply.started":"2022-08-11T18:27:08.451142Z","shell.execute_reply":"2022-08-11T18:27:08.461340Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_exout = df_exout[(df_exout.SaleCondition != 'Abnorml')]\ndf_exout = df_exout[(df_exout.SaleCondition != 'Family')]\ndf_exout.shape","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:28:05.174980Z","iopub.execute_input":"2022-08-11T18:28:05.175387Z","iopub.status.idle":"2022-08-11T18:28:05.192211Z","shell.execute_reply.started":"2022-08-11T18:28:05.175354Z","shell.execute_reply":"2022-08-11T18:28:05.190902Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Despite all the measures above, 'SalePrice' still very much skewed to the right\nfig, ax = plt.subplots(figsize=(10,5))\nsns.distplot(df_exout.SalePrice, bins=30, kde=True, ax=ax)\nplt.title('Price per Square Foot Histogram', fontsize=15)\nplt.xlabel('$ Price/sq ft', fontsize=12);","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:28:23.194838Z","iopub.execute_input":"2022-08-11T18:28:23.195309Z","iopub.status.idle":"2022-08-11T18:28:23.569600Z","shell.execute_reply.started":"2022-08-11T18:28:23.195271Z","shell.execute_reply":"2022-08-11T18:28:23.568195Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Much less of a right-skew to the price per square foot, but still some\nfig, ax = plt.subplots(figsize=(10,5))\nsns.distplot(df_exout.PriceSF, bins=30, kde=True, ax=ax)\nplt.title('Price per Square Foot Histogram', fontsize=15)\nplt.xlabel('$ Price/sq ft', fontsize=12);","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:29:16.930169Z","iopub.execute_input":"2022-08-11T18:29:16.930762Z","iopub.status.idle":"2022-08-11T18:29:17.345790Z","shell.execute_reply.started":"2022-08-11T18:29:16.930711Z","shell.execute_reply":"2022-08-11T18:29:17.344539Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Log transformation of the target variable","metadata":{}},{"cell_type":"code","source":"# Transform both sale price variables by taking the natural log to reduce the right-skew of the distributions\n# The log transformation is also appropriate given that both variables have only non-zero positive values\ndf_exout['LnSalePrice'] = np.log(df_exout.SalePrice)\ndf_exout['LnPriceSF'] = np.log(df_exout.PriceSF)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:35:25.874234Z","iopub.execute_input":"2022-08-11T18:35:25.874645Z","iopub.status.idle":"2022-08-11T18:35:25.882923Z","shell.execute_reply.started":"2022-08-11T18:35:25.874612Z","shell.execute_reply":"2022-08-11T18:35:25.881564Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# The histogram of the ln SalePrice is clearly more symmetric\nfig, ax = plt.subplots(figsize=(10,5))\nsns.distplot(df_exout.LnSalePrice, bins=50, kde=True, ax=ax, color='forestgreen')\nplt.title('Ln Sales Price Histogram', fontsize=15)\nplt.xlabel('Natural Log of Sales Price ($)', fontsize=12)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:35:50.387101Z","iopub.execute_input":"2022-08-11T18:35:50.387539Z","iopub.status.idle":"2022-08-11T18:35:50.821136Z","shell.execute_reply.started":"2022-08-11T18:35:50.387503Z","shell.execute_reply":"2022-08-11T18:35:50.819915Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Natural log of 'SalePrice' is closer to a normal distribution\nfig, ax = plt.subplots(figsize=(10,6))\nfig = probscale.probplot(df_exout.LnSalePrice, ax=ax, plottype='prob')\n\nax.set_xlim(0.1, 99.9)\nax.set_xscale('prob')\nsns.despine(fig=fig)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:36:13.298841Z","iopub.execute_input":"2022-08-11T18:36:13.299279Z","iopub.status.idle":"2022-08-11T18:36:13.839428Z","shell.execute_reply.started":"2022-08-11T18:36:13.299244Z","shell.execute_reply":"2022-08-11T18:36:13.838202Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"sns.pairplot(df_exout[['SalePrice', 'PriceSF', 'GrLivArea', 'BaseLivArea', 'Age']])","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:36:27.116904Z","iopub.execute_input":"2022-08-11T18:36:27.117324Z","iopub.status.idle":"2022-08-11T18:36:32.022446Z","shell.execute_reply.started":"2022-08-11T18:36:27.117289Z","shell.execute_reply":"2022-08-11T18:36:32.020944Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Creating a \"location\" variable","metadata":{}},{"cell_type":"code","source":"# Examining the 'PriceSF' for the various neighbourhoods, it is obvious that the price premium differs significantly\n# across them\nneigh_mean = df_exout['PriceSF'].groupby(df_exout['Neighborhood']).count().sort_values()\nneigh_mean","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:37:28.329558Z","iopub.execute_input":"2022-08-11T18:37:28.330059Z","iopub.status.idle":"2022-08-11T18:37:28.343807Z","shell.execute_reply.started":"2022-08-11T18:37:28.330016Z","shell.execute_reply":"2022-08-11T18:37:28.342618Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"len(neigh_mean)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:37:43.229540Z","iopub.execute_input":"2022-08-11T18:37:43.230346Z","iopub.status.idle":"2022-08-11T18:37:43.238036Z","shell.execute_reply.started":"2022-08-11T18:37:43.230305Z","shell.execute_reply":"2022-08-11T18:37:43.236632Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Certain clusters of neighbourhoods do share similar distributions of 'SalePrice' and/or 'PriceSF'\nfig, ax = plt.subplots(figsize=(15,6))\n\nsns.stripplot(x = df_exout.Neighborhood, y = df_exout.SalePrice, order = np.sort(df_exout.Neighborhood.unique()),\n              jitter=0.1, alpha=0.5, ax=ax)\nplt.title('Neighborhood Stripplot', fontsize=24)\nplt.xticks(rotation=45)\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:38:04.147395Z","iopub.execute_input":"2022-08-11T18:38:04.147780Z","iopub.status.idle":"2022-08-11T18:38:04.847063Z","shell.execute_reply.started":"2022-08-11T18:38:04.147749Z","shell.execute_reply":"2022-08-11T18:38:04.845831Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Creating a new numeric ordinal variable for 'Functional'\ndef functional_numeric(x):\n    if 'Typ' in x:\n        return 8\n    elif 'Min1' in x:\n        return 7\n    elif 'Min2' in x:\n        return 6\n    elif 'Mod' in x:\n        return 5\n    elif 'Maj1' in x:\n        return 4\n    elif 'Maj2' in x:\n        return 3\n    elif 'Sev' in x:\n        return 2    \n    else:\n        return 1\n    \ndf_exout['Functional_Num'] = df_exout.Functional.map(functional_numeric)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:38:48.415014Z","iopub.execute_input":"2022-08-11T18:38:48.415428Z","iopub.status.idle":"2022-08-11T18:38:48.424991Z","shell.execute_reply.started":"2022-08-11T18:38:48.415395Z","shell.execute_reply":"2022-08-11T18:38:48.423582Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Creating a new numeric ordinal variable for external conditon or 'ExterCond'\ndef extercond_numeric(x):\n    if 'Ex' in x:\n        return 5\n    elif 'Gd' in x:\n        return 4\n    elif 'TA' in x:\n        return 3\n    elif 'Fa' in x:\n        return 2\n    else:\n        return 1\n    \ndf_exout['ExterCond_Num'] = df_exout.ExterCond.map(extercond_numeric)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:39:02.633003Z","iopub.execute_input":"2022-08-11T18:39:02.633845Z","iopub.status.idle":"2022-08-11T18:39:02.645353Z","shell.execute_reply.started":"2022-08-11T18:39:02.633794Z","shell.execute_reply":"2022-08-11T18:39:02.644162Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Creating a new numeric ordinal variable for external quality or 'ExterQual'\ndef exterqual_numeric(x):\n    if 'Ex' in x:\n        return 5\n    elif 'Gd' in x:\n        return 4\n    elif 'TA' in x:\n        return 3\n    elif 'Fa' in x:\n        return 2\n    else:\n        return 1\n    \ndf_exout['ExterQual_Num'] = df_exout.ExterQual.map(exterqual_numeric)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:39:21.454788Z","iopub.execute_input":"2022-08-11T18:39:21.455282Z","iopub.status.idle":"2022-08-11T18:39:21.465283Z","shell.execute_reply.started":"2022-08-11T18:39:21.455240Z","shell.execute_reply":"2022-08-11T18:39:21.463924Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"print(df_exout['OverallQual'].mean())\nprint(df_exout['OverallCond'].mean())\nprint(df_exout['ExterQual_Num'].mean())\nprint(df_exout['ExterCond_Num'].mean())\nprint(df_exout['Functional_Num'].mean())","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:39:34.074892Z","iopub.execute_input":"2022-08-11T18:39:34.075304Z","iopub.status.idle":"2022-08-11T18:39:34.086488Z","shell.execute_reply.started":"2022-08-11T18:39:34.075264Z","shell.execute_reply":"2022-08-11T18:39:34.085130Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Construct the location desirability proxy from the 5 individual mean-standardised variables listed above. The construction methodology implies a stronger weighting to the external quality/condition of the building.","metadata":{}},{"cell_type":"code","source":"df_exout['Location'] = ((df_exout['OverallQual']/df_exout['OverallQual'].mean()) \n                        + (df_exout['OverallCond']/df_exout['OverallCond'].mean())\n                        + (df_exout['ExterQual_Num']/df_exout['ExterQual_Num'].mean())\n                        + (df_exout['ExterCond_Num']/df_exout['ExterCond_Num'].mean()) \n                        + (df_exout['Functional_Num']/df_exout['Functional_Num'].mean()))","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:39:58.843989Z","iopub.execute_input":"2022-08-11T18:39:58.844402Z","iopub.status.idle":"2022-08-11T18:39:58.854566Z","shell.execute_reply.started":"2022-08-11T18:39:58.844368Z","shell.execute_reply":"2022-08-11T18:39:58.853762Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_exout.Location.describe()","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:40:09.071884Z","iopub.execute_input":"2022-08-11T18:40:09.072310Z","iopub.status.idle":"2022-08-11T18:40:09.083178Z","shell.execute_reply.started":"2022-08-11T18:40:09.072271Z","shell.execute_reply":"2022-08-11T18:40:09.082291Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_exout.Location.median()","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:40:23.532740Z","iopub.execute_input":"2022-08-11T18:40:23.533187Z","iopub.status.idle":"2022-08-11T18:40:23.541665Z","shell.execute_reply.started":"2022-08-11T18:40:23.533153Z","shell.execute_reply":"2022-08-11T18:40:23.540320Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# The average 'Location' scores of the respective neighbourhoods\ndf_exout['Location'].groupby([df_exout.Neighborhood]).mean().sort_values()","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:40:38.135563Z","iopub.execute_input":"2022-08-11T18:40:38.136123Z","iopub.status.idle":"2022-08-11T18:40:38.150344Z","shell.execute_reply.started":"2022-08-11T18:40:38.136074Z","shell.execute_reply":"2022-08-11T18:40:38.149249Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_exout['SalePrice'].groupby([df_exout.Neighborhood]).count().sort_values()","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:40:53.000920Z","iopub.execute_input":"2022-08-11T18:40:53.001309Z","iopub.status.idle":"2022-08-11T18:40:53.012685Z","shell.execute_reply.started":"2022-08-11T18:40:53.001276Z","shell.execute_reply":"2022-08-11T18:40:53.011782Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def add_location(x):\n    if 'MeadowV' in x or 'Edwards' in x or 'Sawyer' in x or 'Landmrk' in x or 'SWISU' in x or 'BrDale' in x or 'IDOTRR' in x:\n        return 1\n    elif 'NAmes' in x or 'Mitchel' in x or 'BrkSide' in x or 'NPkVill' in x or 'OldTown' in x or 'ClearCr' in x or 'Gilbert' in x:\n        return 2\n    elif 'SawyerW' in x or 'NWAmes' in x or 'Crawfor' in x or 'CollgCr' in x or 'Blueste' in x or 'GrnHill' in x or 'Blmngtn' in x:\n        return 3\n    else:\n        return 4","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:41:54.676172Z","iopub.execute_input":"2022-08-11T18:41:54.676654Z","iopub.status.idle":"2022-08-11T18:41:54.685295Z","shell.execute_reply.started":"2022-08-11T18:41:54.676617Z","shell.execute_reply":"2022-08-11T18:41:54.684006Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_exout['Location'] = df_exout.Neighborhood.map(add_location)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:42:06.768102Z","iopub.execute_input":"2022-08-11T18:42:06.768516Z","iopub.status.idle":"2022-08-11T18:42:06.776830Z","shell.execute_reply.started":"2022-08-11T18:42:06.768484Z","shell.execute_reply":"2022-08-11T18:42:06.775568Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_exout['Location'].value_counts()","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:42:14.625962Z","iopub.execute_input":"2022-08-11T18:42:14.626366Z","iopub.status.idle":"2022-08-11T18:42:14.635676Z","shell.execute_reply.started":"2022-08-11T18:42:14.626333Z","shell.execute_reply":"2022-08-11T18:42:14.634446Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Positive correlation between 'SalePrice' and 'PriceSF' with 'Location'\nprint(df_exout['SalePrice'].groupby(df_exout.Location).mean())\nprint(df_exout['PriceSF'].groupby(df_exout.Location).mean())","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:42:26.374799Z","iopub.execute_input":"2022-08-11T18:42:26.375213Z","iopub.status.idle":"2022-08-11T18:42:26.388446Z","shell.execute_reply.started":"2022-08-11T18:42:26.375181Z","shell.execute_reply":"2022-08-11T18:42:26.387275Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Total mean house size also rises with 'Location'\nprint(df_exout['GrLivArea'].groupby(df_exout.Location).mean())\nprint(df_exout['GrLivArea'].groupby(df_exout.Location).median())","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:42:54.851756Z","iopub.execute_input":"2022-08-11T18:42:54.852938Z","iopub.status.idle":"2022-08-11T18:42:54.865554Z","shell.execute_reply.started":"2022-08-11T18:42:54.852894Z","shell.execute_reply":"2022-08-11T18:42:54.864564Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"So it does appear that the 'Location' score does tend to be positively correlated with house prices, and somewhat with property sizes too.","metadata":{}},{"cell_type":"markdown","source":"### Other feature engineering on sale year, seasonality, zoning, and proximinity to railways & artery roads","metadata":{}},{"cell_type":"code","source":"df_exout['YrSold'] = df_exout['YrSold'].astype(np.int64)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:43:31.187925Z","iopub.execute_input":"2022-08-11T18:43:31.188378Z","iopub.status.idle":"2022-08-11T18:43:31.195707Z","shell.execute_reply.started":"2022-08-11T18:43:31.188344Z","shell.execute_reply":"2022-08-11T18:43:31.194449Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Looking at the 'SalePrice' and 'PriceSF' grouped by 'YrSold', there is an obvious impact of the financial crisis \n# that hit in 2007-2008. The mean 'SoldPrice' and 'PriceSF' in 2009 and 2010 are still lower than those in 2007.\nprint(df_exout.SalePrice.groupby(df_exout.YrSold).mean())\nprint(df_exout.PriceSF.groupby(df_exout.YrSold).mean())","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:43:45.835593Z","iopub.execute_input":"2022-08-11T18:43:45.836018Z","iopub.status.idle":"2022-08-11T18:43:45.848817Z","shell.execute_reply.started":"2022-08-11T18:43:45.835983Z","shell.execute_reply":"2022-08-11T18:43:45.847937Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# House prices were in an uptrend in 2006-2007, but fell in 2008 and have been moribund since. \ndf_exout.PriceSF.groupby(df_exout.YrSold).mean().pct_change()","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:45:24.096935Z","iopub.execute_input":"2022-08-11T18:45:24.097368Z","iopub.status.idle":"2022-08-11T18:45:24.113668Z","shell.execute_reply.started":"2022-08-11T18:45:24.097333Z","shell.execute_reply":"2022-08-11T18:45:24.112344Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Thus, it might be useful to incorporate year dummies in the regression model\ndf_exout = pd.get_dummies(df_exout, columns=['YrSold'])","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:45:38.543193Z","iopub.execute_input":"2022-08-11T18:45:38.543622Z","iopub.status.idle":"2022-08-11T18:45:38.558150Z","shell.execute_reply.started":"2022-08-11T18:45:38.543586Z","shell.execute_reply":"2022-08-11T18:45:38.557163Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Clear price differentials according to the discrete housing zones\ndf_exout.SalePrice.groupby(df_exout.MSZoning).mean()","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:45:51.042330Z","iopub.execute_input":"2022-08-11T18:45:51.042772Z","iopub.status.idle":"2022-08-11T18:45:51.056982Z","shell.execute_reply.started":"2022-08-11T18:45:51.042733Z","shell.execute_reply":"2022-08-11T18:45:51.055366Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Creating a new ordinal variable of zoning, corresponding to the mean values from 'MSZoning'\n# The FV observations in 'MSZoning' are a retirement community development, named \"Floating Village\", which appears \n# wealthy retireesto cater to\n\ndef add_zoning(x):\n    if 'RM' in x:\n        return 1\n    elif 'RH' in x:\n        return 2\n    elif 'RL' in x:\n        return 3\n    else:\n        return 4","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:46:04.722181Z","iopub.execute_input":"2022-08-11T18:46:04.722622Z","iopub.status.idle":"2022-08-11T18:46:04.729225Z","shell.execute_reply.started":"2022-08-11T18:46:04.722583Z","shell.execute_reply":"2022-08-11T18:46:04.728028Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_exout['Zoning'] = df_exout.MSZoning.map(add_zoning)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:46:15.788824Z","iopub.execute_input":"2022-08-11T18:46:15.789598Z","iopub.status.idle":"2022-08-11T18:46:15.795903Z","shell.execute_reply.started":"2022-08-11T18:46:15.789560Z","shell.execute_reply":"2022-08-11T18:46:15.794608Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# No strong linear relationship between 'GrLivArea' or 'BaseLivArea' and 'Zoning'\nprint(df_exout.GrLivArea.groupby(df_exout.Zoning).mean())\ndf_exout.BaseLivArea.groupby(df_exout.Zoning).mean()","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:46:24.243075Z","iopub.execute_input":"2022-08-11T18:46:24.243524Z","iopub.status.idle":"2022-08-11T18:46:24.257485Z","shell.execute_reply.started":"2022-08-11T18:46:24.243490Z","shell.execute_reply":"2022-08-11T18:46:24.256249Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# There appears to be a sizeable jump in prices in the low density and FV categories\nprint(df_exout.SalePrice.groupby(df_exout.Zoning).mean())\ndf_exout.PriceSF.groupby(df_exout.Zoning).mean()","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:46:38.096363Z","iopub.execute_input":"2022-08-11T18:46:38.096763Z","iopub.status.idle":"2022-08-11T18:46:38.112531Z","shell.execute_reply.started":"2022-08-11T18:46:38.096730Z","shell.execute_reply":"2022-08-11T18:46:38.111308Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# There is a relationship between 'Location' and 'Zoning' with Locations 3 and 4 \n# primarily in zones 3 and 4 (low density and the retirement village)\ndf_exout.Location.groupby(df_exout.Zoning).value_counts()","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:46:58.876587Z","iopub.execute_input":"2022-08-11T18:46:58.877090Z","iopub.status.idle":"2022-08-11T18:46:58.890690Z","shell.execute_reply.started":"2022-08-11T18:46:58.877050Z","shell.execute_reply":"2022-08-11T18:46:58.889577Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# There is a general increase inprices as we go from Zone 1 to 3, especially in mean 'PriceSF'\nfig, ax = plt.subplots(ncols=2, figsize=(14, 6))\n\nax[0].scatter(df_exout.SalePrice, df_exout.Zoning)\nax[0].set_title('Sale Price, $', fontsize=16)\nax[1].scatter(df_exout.PriceSF, df_exout.Zoning)\nax[1].set_title('Price per Square Feet, $', fontsize=16)\n#fig.suptitle('Sale Prices by Zone', fontsize=24)\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:47:11.216215Z","iopub.execute_input":"2022-08-11T18:47:11.216655Z","iopub.status.idle":"2022-08-11T18:47:11.689158Z","shell.execute_reply.started":"2022-08-11T18:47:11.216617Z","shell.execute_reply":"2022-08-11T18:47:11.688052Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Analysis above suggests that it might be more useful to model 'Zoning' using dummies instead of an ordinal variable\ndf_exout['Zone_ordinal'] = df_exout['Zoning']\ndf_exout = pd.get_dummies(df_exout, columns=['Zoning'])","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:47:57.255367Z","iopub.execute_input":"2022-08-11T18:47:57.255805Z","iopub.status.idle":"2022-08-11T18:47:57.273397Z","shell.execute_reply.started":"2022-08-11T18:47:57.255762Z","shell.execute_reply":"2022-08-11T18:47:57.271945Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# The overwhelming majority of observations fall into the double-'Norm' category in terms of\n# Condition1 and Condition2, which is a measure of environmental condition\ndf_exout.Condition1.groupby(df_exout.Condition2).value_counts()","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:48:09.327963Z","iopub.execute_input":"2022-08-11T18:48:09.328399Z","iopub.status.idle":"2022-08-11T18:48:09.340804Z","shell.execute_reply.started":"2022-08-11T18:48:09.328362Z","shell.execute_reply":"2022-08-11T18:48:09.339935Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Previous studies have shown a negative relationship between housing prices and being adjacent \n# or near to a major roadway or a railway line.\n# Creating variables to capture this negative enviromental impact\ndef add_roadrail1(x):\n    if 'Artery' in x:\n        return 1\n    elif 'RRAn' in x:\n        return 1\n    elif 'RRNn' in x:\n        return 1\n    elif 'RRAe' in x:\n        return 1\n    elif 'RRNe' in x:\n        return 1\n    else:\n        return 0\n\ndf_exout['RoadRail1'] = df_exout.Condition1.map(add_roadrail1)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:49:28.831462Z","iopub.execute_input":"2022-08-11T18:49:28.831975Z","iopub.status.idle":"2022-08-11T18:49:28.842197Z","shell.execute_reply.started":"2022-08-11T18:49:28.831935Z","shell.execute_reply":"2022-08-11T18:49:28.840895Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def add_roadrail2(x):\n    if 'Artery' in x:\n        return 1\n    elif 'RRAn' in x:\n        return 1\n    elif 'RRNn' in x:\n        return 1\n    elif 'RRAe' in x:\n        return 1\n    elif 'RRNe' in x:\n        return 1\n    else:\n        return 0\n\ndf_exout['RoadRail2'] = df_exout.Condition2.map(add_roadrail2)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:49:44.050753Z","iopub.execute_input":"2022-08-11T18:49:44.051190Z","iopub.status.idle":"2022-08-11T18:49:44.062104Z","shell.execute_reply.started":"2022-08-11T18:49:44.051155Z","shell.execute_reply":"2022-08-11T18:49:44.061127Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_exout.RoadRail1.groupby(df_exout.RoadRail2).value_counts()","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:49:58.202768Z","iopub.execute_input":"2022-08-11T18:49:58.203215Z","iopub.status.idle":"2022-08-11T18:49:58.219235Z","shell.execute_reply.started":"2022-08-11T18:49:58.203179Z","shell.execute_reply":"2022-08-11T18:49:58.217783Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Combining the 2 proximity to artery roads and railroad dummy variables into one\ndf_exout['RoadRail'] = df_exout.apply(lambda x: 1 if (x['RoadRail1'] == 1 | x['RoadRail2'] == 1) \n                                      else 0, axis=1)\ndf_exout.drop(['RoadRail1', 'RoadRail2'], axis=1, inplace=True)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:50:12.114372Z","iopub.execute_input":"2022-08-11T18:50:12.114778Z","iopub.status.idle":"2022-08-11T18:50:12.155725Z","shell.execute_reply.started":"2022-08-11T18:50:12.114744Z","shell.execute_reply":"2022-08-11T18:50:12.154104Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# The RoadRail variable seems to have an impact on mean 'SalePrice' and 'PriceSF' in every 'Location'\nprint(df_exout.SalePrice.groupby([df_exout.RoadRail, df_exout.Location]).mean())\ndf_exout.PriceSF.groupby([df_exout.RoadRail, df_exout.Location]).mean()","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:50:21.932052Z","iopub.execute_input":"2022-08-11T18:50:21.932445Z","iopub.status.idle":"2022-08-11T18:50:21.954018Z","shell.execute_reply.started":"2022-08-11T18:50:21.932415Z","shell.execute_reply":"2022-08-11T18:50:21.952472Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### More (minor) feature engineering","metadata":{}},{"cell_type":"code","source":"# Transforming the 'CentralAir' discrete variable to numeric\ndf_exout['CentralAirNum'] = df_exout.apply(lambda x: 1 if (x['CentralAir'] == 'Y') \n                                                           else 0, axis=1)\n#df_exout[['CentralAir', 'CentralAirNum']]\ndf_exout['CentralAirNum'].value_counts()","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:50:46.080944Z","iopub.execute_input":"2022-08-11T18:50:46.081679Z","iopub.status.idle":"2022-08-11T18:50:46.113272Z","shell.execute_reply.started":"2022-08-11T18:50:46.081629Z","shell.execute_reply":"2022-08-11T18:50:46.112024Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Creating variables on positive neighbourhood amenities\ndef add_amenities1(x):\n    if 'PosN' in x:\n        return 1\n    elif 'PosA' in x:\n        return 1\n    else:\n        return 0\n\ndf_exout['Amenities1'] = df_exout.Condition1.map(add_amenities1)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:50:57.927052Z","iopub.execute_input":"2022-08-11T18:50:57.928041Z","iopub.status.idle":"2022-08-11T18:50:57.936257Z","shell.execute_reply.started":"2022-08-11T18:50:57.927999Z","shell.execute_reply":"2022-08-11T18:50:57.935278Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def add_amenities2(x):\n    if 'PosN' in x:\n        return 1\n    elif 'PosA' in x:\n        return 1\n    else:\n        return 0\n\ndf_exout['Amenities2'] = df_exout.Condition2.map(add_amenities2)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:51:08.749296Z","iopub.execute_input":"2022-08-11T18:51:08.749761Z","iopub.status.idle":"2022-08-11T18:51:08.759220Z","shell.execute_reply.started":"2022-08-11T18:51:08.749701Z","shell.execute_reply":"2022-08-11T18:51:08.758038Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Combining the amenities dummy variables into one\ndf_exout['Amenities'] = df_exout.apply(lambda x: 1 if (x['Amenities1'] == 1 | x['Amenities2'] == 1) \n                                      else 0, axis=1)\ndf_exout.drop(['Amenities1', 'Amenities2'], axis=1, inplace=True)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:51:18.551396Z","iopub.execute_input":"2022-08-11T18:51:18.551803Z","iopub.status.idle":"2022-08-11T18:51:18.595934Z","shell.execute_reply.started":"2022-08-11T18:51:18.551769Z","shell.execute_reply":"2022-08-11T18:51:18.594634Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# 'Amenities'=1 does seem to raise the mean 'SalePrice' but has mixed impact on 'PriceSF'\nprint(df_exout.SalePrice.groupby([df_exout.Amenities, df_exout.Location]).mean())\ndf_exout.PriceSF.groupby([df_exout.Amenities, df_exout.Location]).mean()","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:51:29.330224Z","iopub.execute_input":"2022-08-11T18:51:29.331134Z","iopub.status.idle":"2022-08-11T18:51:29.351674Z","shell.execute_reply.started":"2022-08-11T18:51:29.331093Z","shell.execute_reply":"2022-08-11T18:51:29.350349Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# At this point, all the observations in the sample have full public utility service\ndf_exout[(df_exout.Utilities != 'AllPub')].count()","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:51:43.826983Z","iopub.execute_input":"2022-08-11T18:51:43.827435Z","iopub.status.idle":"2022-08-11T18:51:43.838711Z","shell.execute_reply.started":"2022-08-11T18:51:43.827399Z","shell.execute_reply":"2022-08-11T18:51:43.837518Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Adding up the above ground bathrooms and assigning to a new variable\ndf_exout['Bathrooms'] = df_exout.FullBath + (0.5*df_exout.HalfBath)\ndf_exout.Bathrooms.describe()","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:52:05.521810Z","iopub.execute_input":"2022-08-11T18:52:05.522263Z","iopub.status.idle":"2022-08-11T18:52:05.536126Z","shell.execute_reply.started":"2022-08-11T18:52:05.522228Z","shell.execute_reply":"2022-08-11T18:52:05.535078Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Positive relationship between the overall quality of the housing structure and the 'SalePrice' and 'PriceSF'\nprint(df_exout.SalePrice.groupby([df_exout.OverallQual]).mean())\ndf_exout.PriceSF.groupby([df_exout.OverallQual]).mean()","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:52:15.642258Z","iopub.execute_input":"2022-08-11T18:52:15.642732Z","iopub.status.idle":"2022-08-11T18:52:15.660522Z","shell.execute_reply.started":"2022-08-11T18:52:15.642691Z","shell.execute_reply":"2022-08-11T18:52:15.659336Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Interestingly 'OverallQual' is not well correlated with 'OverallCond'\ndf_exout.OverallQual.corr(df_exout.OverallCond)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:52:28.558545Z","iopub.execute_input":"2022-08-11T18:52:28.558961Z","iopub.status.idle":"2022-08-11T18:52:28.566682Z","shell.execute_reply.started":"2022-08-11T18:52:28.558926Z","shell.execute_reply":"2022-08-11T18:52:28.565900Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# But it is moderately correlated with the 'Location' variable as one would expect\ndf_exout.OverallQual.corr(df_exout.Location)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:52:37.149053Z","iopub.execute_input":"2022-08-11T18:52:37.149454Z","iopub.status.idle":"2022-08-11T18:52:37.157703Z","shell.execute_reply.started":"2022-08-11T18:52:37.149421Z","shell.execute_reply":"2022-08-11T18:52:37.156885Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Overall Condition of the house appears to have mixed impact on prices\nprint(df_exout.SalePrice.groupby([df_exout.OverallCond]).mean())\ndf_exout.PriceSF.groupby([df_exout.OverallCond]).mean()","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:52:47.667147Z","iopub.execute_input":"2022-08-11T18:52:47.668383Z","iopub.status.idle":"2022-08-11T18:52:47.685304Z","shell.execute_reply.started":"2022-08-11T18:52:47.668332Z","shell.execute_reply":"2022-08-11T18:52:47.684373Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Dummifying other potential important features","metadata":{}},{"cell_type":"code","source":"df_exout = pd.get_dummies(df_exout, columns=['Street'])","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:53:59.398277Z","iopub.execute_input":"2022-08-11T18:53:59.398756Z","iopub.status.idle":"2022-08-11T18:53:59.415147Z","shell.execute_reply.started":"2022-08-11T18:53:59.398721Z","shell.execute_reply":"2022-08-11T18:53:59.413618Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_exout.drop(['Street_Grvl'], axis=1, inplace=True)\ndf_exout.rename(columns={'Street_Pave': 'Street_dum'}, inplace=True)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:54:06.421414Z","iopub.execute_input":"2022-08-11T18:54:06.422099Z","iopub.status.idle":"2022-08-11T18:54:06.433056Z","shell.execute_reply.started":"2022-08-11T18:54:06.422060Z","shell.execute_reply":"2022-08-11T18:54:06.431969Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_exout = pd.get_dummies(df_exout, columns=['KitchenQual'])","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:54:18.463804Z","iopub.execute_input":"2022-08-11T18:54:18.464261Z","iopub.status.idle":"2022-08-11T18:54:18.481701Z","shell.execute_reply.started":"2022-08-11T18:54:18.464224Z","shell.execute_reply":"2022-08-11T18:54:18.480608Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Creating a dummy variable for flat roofs\ndf_exout['FlatRoof_dum'] = df_exout.RoofStyle.apply(lambda x: 1 if x=='Flat' else 0)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:54:26.931016Z","iopub.execute_input":"2022-08-11T18:54:26.931716Z","iopub.status.idle":"2022-08-11T18:54:26.939527Z","shell.execute_reply.started":"2022-08-11T18:54:26.931666Z","shell.execute_reply":"2022-08-11T18:54:26.938536Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Creating a dummy variable for garage\ndf_exout['Garage_dum'] = df_exout.GarageQual.apply(lambda x: 0 if pd.isnull(x)==True else 1)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:54:39.569091Z","iopub.execute_input":"2022-08-11T18:54:39.569476Z","iopub.status.idle":"2022-08-11T18:54:39.578700Z","shell.execute_reply.started":"2022-08-11T18:54:39.569445Z","shell.execute_reply":"2022-08-11T18:54:39.577555Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Creating a dummy variable for flat property contour\ndf_exout['FlatContour_dum'] = df_exout.LandContour.apply(lambda x: 1 if x=='Lvl' else 0)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:54:49.595816Z","iopub.execute_input":"2022-08-11T18:54:49.596318Z","iopub.status.idle":"2022-08-11T18:54:49.604626Z","shell.execute_reply.started":"2022-08-11T18:54:49.596269Z","shell.execute_reply":"2022-08-11T18:54:49.603225Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Creating a dummy variable for houses higher than one storey\ndf_exout['TwoStory_dum'] = df_exout['2ndFlrSF'].apply(lambda x: 1 if x>0 else 0)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:54:58.809968Z","iopub.execute_input":"2022-08-11T18:54:58.810442Z","iopub.status.idle":"2022-08-11T18:54:58.818697Z","shell.execute_reply.started":"2022-08-11T18:54:58.810404Z","shell.execute_reply":"2022-08-11T18:54:58.817604Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_exout.shape","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:55:07.266811Z","iopub.execute_input":"2022-08-11T18:55:07.267276Z","iopub.status.idle":"2022-08-11T18:55:07.275375Z","shell.execute_reply.started":"2022-08-11T18:55:07.267240Z","shell.execute_reply":"2022-08-11T18:55:07.274113Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Finalising the dataset","metadata":{}},{"cell_type":"code","source":"pd.options.display.max_columns = None\ndf_exout.head()","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:55:31.668958Z","iopub.execute_input":"2022-08-11T18:55:31.669393Z","iopub.status.idle":"2022-08-11T18:55:31.756256Z","shell.execute_reply.started":"2022-08-11T18:55:31.669357Z","shell.execute_reply":"2022-08-11T18:55:31.755170Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_exout = df_exout[['SalePrice', 'LnSalePrice', 'Age', 'GrLivArea', 'BaseLivArea', 'Location', 'Amenities', \n                     'RoadRail', 'BedroomAbvGr', 'Bathrooms', 'OverallCond', 'OverallQual', 'LotFrontage', \n                     'LotArea', 'TwoStory_dum', 'FlatContour_dum', 'FlatRoof_dum', 'GarageArea', 'Garage_dum', \n                     'CentralAirNum', 'LowQualFinSF', 'Fireplaces', 'KitchenQual_Ex', 'Zoning_2', 'Zoning_3', \n                     'Zoning_4', 'YrSold_2007', 'YrSold_2008', 'YrSold_2009', 'YrSold_2010']]","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:55:52.064546Z","iopub.execute_input":"2022-08-11T18:55:52.064988Z","iopub.status.idle":"2022-08-11T18:55:52.074721Z","shell.execute_reply.started":"2022-08-11T18:55:52.064952Z","shell.execute_reply":"2022-08-11T18:55:52.073384Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Now on to the regression analysis...","metadata":{}},{"cell_type":"code","source":"df = df_exout\ndf.info()","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:58:26.669808Z","iopub.execute_input":"2022-08-11T18:58:26.670201Z","iopub.status.idle":"2022-08-11T18:58:26.685206Z","shell.execute_reply.started":"2022-08-11T18:58:26.670170Z","shell.execute_reply":"2022-08-11T18:58:26.684219Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Deleting the single null value in 'GarageArea'\ndf = df[~df['GarageArea'].isnull()]","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:58:59.805329Z","iopub.execute_input":"2022-08-11T18:58:59.805760Z","iopub.status.idle":"2022-08-11T18:58:59.813155Z","shell.execute_reply.started":"2022-08-11T18:58:59.805718Z","shell.execute_reply":"2022-08-11T18:58:59.812040Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df.isnull().values.any()","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:59:07.829878Z","iopub.execute_input":"2022-08-11T18:59:07.830723Z","iopub.status.idle":"2022-08-11T18:59:07.836728Z","shell.execute_reply.started":"2022-08-11T18:59:07.830688Z","shell.execute_reply":"2022-08-11T18:59:07.836007Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"fig, ax = plt.subplots(figsize=(18, 16))\n\nsns.set(font_scale=0.8)\nsns.heatmap(df.corr(), annot=True, cmap='coolwarm', ax=ax)\nax.set_title(\"Correlation Matrix of Ames Variables\", fontsize=16)\nax.set_xticklabels(ax.get_xmajorticklabels(), fontsize=8)\nax.set_yticklabels(ax.get_ymajorticklabels(), fontsize=8)\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2022-08-11T18:59:19.750602Z","iopub.execute_input":"2022-08-11T18:59:19.751083Z","iopub.status.idle":"2022-08-11T18:59:24.179389Z","shell.execute_reply.started":"2022-08-11T18:59:19.751043Z","shell.execute_reply":"2022-08-11T18:59:24.178278Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"The correlation matrix shows that the explanatory variables tend not to have high correlations with each other. The explanatory variable pairs that exhibit correlations is bigger than 0.70 is GrLivArea-Bathrooms. Bathrooms also has elevated correlation with Age. Given this, I decided to drop Bathrooms from the regression analysis.","metadata":{}},{"cell_type":"code","source":"df = df.drop(['Bathrooms'], axis=1)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:00:43.941962Z","iopub.execute_input":"2022-08-11T19:00:43.942400Z","iopub.status.idle":"2022-08-11T19:00:43.949767Z","shell.execute_reply.started":"2022-08-11T19:00:43.942361Z","shell.execute_reply":"2022-08-11T19:00:43.948518Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Dividing the data into pre-2010 and a 2010 holdout","metadata":{}},{"cell_type":"code","source":"# Assigning the 2006-2009 data to another dataset\ndf_0609 = df.loc[df['YrSold_2010'] != 1]\ndf_0609.shape","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:01:11.440249Z","iopub.execute_input":"2022-08-11T19:01:11.440712Z","iopub.status.idle":"2022-08-11T19:01:11.451345Z","shell.execute_reply.started":"2022-08-11T19:01:11.440668Z","shell.execute_reply":"2022-08-11T19:01:11.450144Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Assigning 2010 data as the holdout test set\ndf_2010 = df.loc[df['YrSold_2010'] == 1]","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:01:20.170202Z","iopub.execute_input":"2022-08-11T19:01:20.170646Z","iopub.status.idle":"2022-08-11T19:01:20.177930Z","shell.execute_reply.started":"2022-08-11T19:01:20.170609Z","shell.execute_reply":"2022-08-11T19:01:20.176648Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Functions","metadata":{}},{"cell_type":"code","source":"# Function for scoring training set\ndef train_scores(model, X, y):\n    '''\n    model: fitted model\n    X: Matrix of explanatory variables (train set)\n    y: Dependant variable (train set)\n    '''\n    cv_scores = cross_val_score(model, X, y, cv=5) # 5-fold cross-validation\n\n    print('Training Score:', np.round(model.score(X, y), 4))\n    print('Cross-validation scores:', np.round(cv_scores, 4))\n    print('Mean cross-validation score:', np.round(cv_scores.mean(), 4))","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:05:18.478954Z","iopub.execute_input":"2022-08-11T19:05:18.480080Z","iopub.status.idle":"2022-08-11T19:05:18.487344Z","shell.execute_reply.started":"2022-08-11T19:05:18.480033Z","shell.execute_reply":"2022-08-11T19:05:18.486053Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Function for scoring test set\ndef test_scores(model, X, y):\n    '''\n    model: fitted model\n    X: Matrix of explanatory variables (test set)\n    y: Dependant variable (test set)\n    '''\n    print('Test Score:', np.round(model.score(X, y), 4))","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:05:34.931233Z","iopub.execute_input":"2022-08-11T19:05:34.931630Z","iopub.status.idle":"2022-08-11T19:05:34.937139Z","shell.execute_reply.started":"2022-08-11T19:05:34.931598Z","shell.execute_reply":"2022-08-11T19:05:34.935880Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Function for MSE & RMSE scoring\ndef accuracy_scores(model, X, y):\n    '''\n    model: fitted model\n    X: Matrix of explanatory variables (test set)\n    y: Dependant variable (test set)\n    '''\n    yhat = model.predict(X)\n    print('Mean Squared Error:', np.round(metrics.mean_squared_error(y, yhat), 4)) \n    print('Root Mean Squared Error:', np.round((metrics.mean_squared_error(y, yhat))**0.5, 4)) ","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:05:45.414521Z","iopub.execute_input":"2022-08-11T19:05:45.414929Z","iopub.status.idle":"2022-08-11T19:05:45.422547Z","shell.execute_reply.started":"2022-08-11T19:05:45.414894Z","shell.execute_reply":"2022-08-11T19:05:45.421202Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Function for plotting histogram of residuals\ndef resid_histogram(model, X, y, period=''):\n    '''\n    model: fitted model\n    X: Matrix of explanatory variables\n    y: Dependant variable\n    period: String describing data coverage period\n    '''\n\n    yhat = model.predict(X)\n    residuals = y - yhat\n\n    fig, ax = plt.subplots(figsize=(10,6))\n    sns.distplot(residuals, bins=50, kde=True, ax=ax)\n    plt.title(f'OLS Residuals, {period}', fontsize=18)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:05:59.562050Z","iopub.execute_input":"2022-08-11T19:05:59.562521Z","iopub.status.idle":"2022-08-11T19:05:59.569779Z","shell.execute_reply.started":"2022-08-11T19:05:59.562481Z","shell.execute_reply":"2022-08-11T19:05:59.568633Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Regression using Ln SalePrice target on pre-2010","metadata":{}},{"cell_type":"code","source":"y_SP = df_0609['SalePrice']\ny_lnSP = df_0609['LnSalePrice']","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:06:31.571765Z","iopub.execute_input":"2022-08-11T19:06:31.572207Z","iopub.status.idle":"2022-08-11T19:06:31.579168Z","shell.execute_reply.started":"2022-08-11T19:06:31.572172Z","shell.execute_reply":"2022-08-11T19:06:31.577739Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"X = df_0609.drop(['SalePrice', 'LnSalePrice'], axis=1)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:06:43.593992Z","iopub.execute_input":"2022-08-11T19:06:43.594408Z","iopub.status.idle":"2022-08-11T19:06:43.600671Z","shell.execute_reply.started":"2022-08-11T19:06:43.594373Z","shell.execute_reply":"2022-08-11T19:06:43.599869Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"X.shape","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:06:51.700554Z","iopub.execute_input":"2022-08-11T19:06:51.700984Z","iopub.status.idle":"2022-08-11T19:06:51.708251Z","shell.execute_reply.started":"2022-08-11T19:06:51.700946Z","shell.execute_reply":"2022-08-11T19:06:51.707139Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Train-test split of the 2006-2009 data\nX_train, X_test, y_train, y_test = train_test_split(X, y_lnSP, test_size=0.3, random_state=8)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:07:03.034732Z","iopub.execute_input":"2022-08-11T19:07:03.035162Z","iopub.status.idle":"2022-08-11T19:07:03.043817Z","shell.execute_reply.started":"2022-08-11T19:07:03.035130Z","shell.execute_reply":"2022-08-11T19:07:03.042700Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"scaler = StandardScaler()\nX_train = pd.DataFrame(scaler.fit_transform(X_train), columns=X_train.columns)\nX_test = pd.DataFrame(scaler.transform(X_test), columns=X_test.columns)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:07:18.803016Z","iopub.execute_input":"2022-08-11T19:07:18.803439Z","iopub.status.idle":"2022-08-11T19:07:18.819008Z","shell.execute_reply.started":"2022-08-11T19:07:18.803405Z","shell.execute_reply":"2022-08-11T19:07:18.817770Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Fitting ordinary linear regression & getting parameter estimates\nols = LinearRegression()\nols.fit(X_train, y_train)\n\nprint(\"Intercept:\", ols.intercept_)\nprint(\"Coefficients:\", ols.coef_)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:07:30.761337Z","iopub.execute_input":"2022-08-11T19:07:30.761772Z","iopub.status.idle":"2022-08-11T19:07:30.789743Z","shell.execute_reply.started":"2022-08-11T19:07:30.761732Z","shell.execute_reply":"2022-08-11T19:07:30.788115Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# OLS training set scores, including CV scores\ntrain_scores(ols, X_train, y_train)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:07:43.254788Z","iopub.execute_input":"2022-08-11T19:07:43.255238Z","iopub.status.idle":"2022-08-11T19:07:43.383976Z","shell.execute_reply.started":"2022-08-11T19:07:43.255202Z","shell.execute_reply":"2022-08-11T19:07:43.382397Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Shuffled 5-fold cross validation scores are rather similar\nkf = KFold(n_splits=5, shuffle=True, random_state=1)\ncv_scores_shuffled = cross_val_score(ols, X_train, y_train, cv=kf)\n\nprint('Shuffled cross validation score:', np.round(cv_scores_shuffled, 4))\nprint('Mean shuffled cross validation score:', np.round(cv_scores_shuffled.mean(), 4))","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:08:05.411062Z","iopub.execute_input":"2022-08-11T19:08:05.411467Z","iopub.status.idle":"2022-08-11T19:08:05.485070Z","shell.execute_reply.started":"2022-08-11T19:08:05.411433Z","shell.execute_reply":"2022-08-11T19:08:05.483351Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# OLS test set score\ntest_scores(ols, X_test, y_test)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:08:21.288092Z","iopub.execute_input":"2022-08-11T19:08:21.288490Z","iopub.status.idle":"2022-08-11T19:08:21.307879Z","shell.execute_reply.started":"2022-08-11T19:08:21.288458Z","shell.execute_reply":"2022-08-11T19:08:21.306034Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# OLS MSE & RMSE scores\naccuracy_scores(ols, X_test, y_test)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:08:34.578472Z","iopub.execute_input":"2022-08-11T19:08:34.578879Z","iopub.status.idle":"2022-08-11T19:08:34.592019Z","shell.execute_reply.started":"2022-08-11T19:08:34.578831Z","shell.execute_reply":"2022-08-11T19:08:34.590494Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Collect the coefficients\ndf_ols_coef = pd.DataFrame(ols.coef_, index=X_train.columns, columns=['Coefficients'])\ndf_ols_coef['Coef_abs'] = df_ols_coef.Coefficients.abs()","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:08:53.551337Z","iopub.execute_input":"2022-08-11T19:08:53.551720Z","iopub.status.idle":"2022-08-11T19:08:53.561456Z","shell.execute_reply.started":"2022-08-11T19:08:53.551688Z","shell.execute_reply":"2022-08-11T19:08:53.559964Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Analysis of the OLS residuals","metadata":{}},{"cell_type":"code","source":"predictions_train = ols.predict(X_train)\npredictions_test = ols.predict(X_test)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:09:47.822518Z","iopub.execute_input":"2022-08-11T19:09:47.823196Z","iopub.status.idle":"2022-08-11T19:09:47.840660Z","shell.execute_reply.started":"2022-08-11T19:09:47.823146Z","shell.execute_reply":"2022-08-11T19:09:47.838950Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Descriptive statistics of training set residuals\nols_residuals_0609 = (y_train - predictions_train)\nols_residuals_0609.describe()","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:09:59.088478Z","iopub.execute_input":"2022-08-11T19:09:59.088893Z","iopub.status.idle":"2022-08-11T19:09:59.101922Z","shell.execute_reply.started":"2022-08-11T19:09:59.088836Z","shell.execute_reply":"2022-08-11T19:09:59.100792Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Acceptable skew and kurtosis values\nprint(\"Skew:\", ols_residuals_0609.skew())\nprint(\"Kurtosis:\", ols_residuals_0609.kurtosis())\nstat, p = shapiro(ols_residuals_0609)\nprint('Shapiro-Wilk test on normality=%.3f, p=%.3f' % (stat, p))","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:10:12.739288Z","iopub.execute_input":"2022-08-11T19:10:12.739723Z","iopub.status.idle":"2022-08-11T19:10:12.747742Z","shell.execute_reply.started":"2022-08-11T19:10:12.739691Z","shell.execute_reply":"2022-08-11T19:10:12.746553Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Histogram of training set residuals show that they are approximately normally distributed with mean 0\n# There is an indication of a left-tail, indicating that the model overpredicts the target variable\n# at the very low end of 'LnSalePrice'\nresid_histogram(ols, X_train, y_train, period='2006-2009 train data')","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:10:24.989546Z","iopub.execute_input":"2022-08-11T19:10:24.990029Z","iopub.status.idle":"2022-08-11T19:10:25.445106Z","shell.execute_reply.started":"2022-08-11T19:10:24.989991Z","shell.execute_reply":"2022-08-11T19:10:25.443869Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"from scipy import stats\nstats.probplot(ols_residuals_0609, dist=\"norm\", plot=plt)\nplt.title(\"Quantile-Quantile Plot, 2006-2009 training set\")","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:10:41.654411Z","iopub.execute_input":"2022-08-11T19:10:41.654846Z","iopub.status.idle":"2022-08-11T19:10:41.941556Z","shell.execute_reply.started":"2022-08-11T19:10:41.654809Z","shell.execute_reply":"2022-08-11T19:10:41.940254Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Descriptive statistics of test set residuals\nols_residuals_test = (y_test - predictions_test)\nols_residuals_test.describe()","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:10:56.815201Z","iopub.execute_input":"2022-08-11T19:10:56.815648Z","iopub.status.idle":"2022-08-11T19:10:56.827746Z","shell.execute_reply.started":"2022-08-11T19:10:56.815610Z","shell.execute_reply":"2022-08-11T19:10:56.826644Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Histogram of test set residuals is symmetric and approximately normal in shape\nresid_histogram(ols, X_test, y_test, period='2006-2009 test data')","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:11:09.758871Z","iopub.execute_input":"2022-08-11T19:11:09.759338Z","iopub.status.idle":"2022-08-11T19:11:10.302419Z","shell.execute_reply.started":"2022-08-11T19:11:09.759302Z","shell.execute_reply":"2022-08-11T19:11:10.301273Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"stats.probplot(ols_residuals_test, dist=\"norm\", plot=plt)\nplt.title(\"Quantile-Quantile Plot, 2006-2009 test set\")","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:11:27.669051Z","iopub.execute_input":"2022-08-11T19:11:27.669489Z","iopub.status.idle":"2022-08-11T19:11:27.966941Z","shell.execute_reply.started":"2022-08-11T19:11:27.669450Z","shell.execute_reply":"2022-08-11T19:11:27.965650Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"ols_residuals_test0609 = (y_test - predictions_test)\nprint(\"Skew:\", ols_residuals_test0609.skew())\nprint(\"Kurtosis:\", ols_residuals_test0609.kurtosis())\nstat, p = shapiro(ols_residuals_test0609)\nprint('Shapiro-Wilk test on normality=%.3f, p=%.3f' % (stat, p))","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:11:44.035156Z","iopub.execute_input":"2022-08-11T19:11:44.035557Z","iopub.status.idle":"2022-08-11T19:11:44.044524Z","shell.execute_reply.started":"2022-08-11T19:11:44.035526Z","shell.execute_reply":"2022-08-11T19:11:44.043100Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Plotting the OLS residuals against the predicted-y and 'GrLivArea'. The residuals appear well-behaved\nfig, ax = plt.subplots(ncols=2, figsize=(15, 6))\nax[0].scatter(ols_residuals_0609, predictions_train)\nax[0].set_title('Residuals vs Predicted Target, 2006-2009 training set', fontsize=14)\nax[1].scatter(ols_residuals_0609, X_train.GrLivArea)\nax[1].set_title('Residuals vs Above Grade Square Footage', fontsize=14)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:12:02.743028Z","iopub.execute_input":"2022-08-11T19:12:02.743424Z","iopub.status.idle":"2022-08-11T19:12:03.321232Z","shell.execute_reply.started":"2022-08-11T19:12:02.743372Z","shell.execute_reply":"2022-08-11T19:12:03.320092Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"fig, ax = plt.subplots(ncols=2, figsize=(15, 6))\nax[0].hist(ols_residuals_test, density=True, bins=30, color='indianred')\nax[0].set_title('OLS Residuals, 2006-2009 test set', fontsize=14)\nax[1].scatter(ols_residuals_test, predictions_test, color='midnightblue')\nax[1].set_title('Residuals vs Predicted Target, 2006-2009 test set', fontsize=14)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:12:26.108069Z","iopub.execute_input":"2022-08-11T19:12:26.108492Z","iopub.status.idle":"2022-08-11T19:12:26.694526Z","shell.execute_reply.started":"2022-08-11T19:12:26.108459Z","shell.execute_reply":"2022-08-11T19:12:26.693278Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Ridge & Lasso regressions","metadata":{}},{"cell_type":"code","source":"# Ridge Cross-Validation\nridge_mod = RidgeCV(alphas=np.logspace(-4, 4, 10), cv=5)\nridge_mod.fit(X_train, y_train)\n\nprint('Best Ridge alpha:', ridge_mod.alpha_)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:13:22.311046Z","iopub.execute_input":"2022-08-11T19:13:22.311806Z","iopub.status.idle":"2022-08-11T19:13:22.861579Z","shell.execute_reply.started":"2022-08-11T19:13:22.311757Z","shell.execute_reply":"2022-08-11T19:13:22.857885Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Ridge training set scores, including CV scores\ntrain_scores(ridge_mod, X_train, y_train)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:13:41.810278Z","iopub.execute_input":"2022-08-11T19:13:41.810711Z","iopub.status.idle":"2022-08-11T19:13:44.791413Z","shell.execute_reply.started":"2022-08-11T19:13:41.810673Z","shell.execute_reply":"2022-08-11T19:13:44.789881Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Ridge test set score\ntest_scores(ridge_mod, X_test, y_test)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:13:57.771240Z","iopub.execute_input":"2022-08-11T19:13:57.771625Z","iopub.status.idle":"2022-08-11T19:13:57.788147Z","shell.execute_reply.started":"2022-08-11T19:13:57.771593Z","shell.execute_reply":"2022-08-11T19:13:57.786470Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Confirmed similar to the above Ridge CV scores\nridge_mod = Ridge(alpha=21.544)\n\nridge_mod.fit(X_train, y_train)\nprint(\"Training Score:\", round(ridge_mod.score(X_train, y_train), 4))\nprint(\"Test Score:\", round(ridge_mod.score(X_test, y_test), 4))","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:14:12.007657Z","iopub.execute_input":"2022-08-11T19:14:12.008290Z","iopub.status.idle":"2022-08-11T19:14:12.060948Z","shell.execute_reply.started":"2022-08-11T19:14:12.008236Z","shell.execute_reply":"2022-08-11T19:14:12.058859Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Ridge MSE & RMSE scores\naccuracy_scores(ridge_mod, X_test, y_test)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:17:31.631840Z","iopub.execute_input":"2022-08-11T19:17:31.633086Z","iopub.status.idle":"2022-08-11T19:17:31.650955Z","shell.execute_reply.started":"2022-08-11T19:17:31.633025Z","shell.execute_reply":"2022-08-11T19:17:31.649074Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Collecting Ridge coefficients\ndf_ridge_coef = pd.DataFrame(ridge_mod.coef_, index=X_train.columns,\n                       columns=['Coefficients'])\ndf_ridge_coef['Coef_abs'] = df_ridge_coef.Coefficients.abs()","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:21:30.019813Z","iopub.execute_input":"2022-08-11T19:21:30.020766Z","iopub.status.idle":"2022-08-11T19:21:30.028159Z","shell.execute_reply.started":"2022-08-11T19:21:30.020713Z","shell.execute_reply":"2022-08-11T19:21:30.027236Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Lasso Cross-Validation\nlasso_mod = LassoCV(alphas=np.logspace(-4, 4, 10), cv=5)\nlasso_mod.fit(X_train, y_train)\n\nprint('Best Lasso alpha:', lasso_mod.alpha_)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:21:46.813987Z","iopub.execute_input":"2022-08-11T19:21:46.814416Z","iopub.status.idle":"2022-08-11T19:21:46.870663Z","shell.execute_reply.started":"2022-08-11T19:21:46.814381Z","shell.execute_reply":"2022-08-11T19:21:46.869087Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Lasso training set scores, including CV scores\ntrain_scores(lasso_mod, X_train, y_train)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:22:00.351104Z","iopub.execute_input":"2022-08-11T19:22:00.351489Z","iopub.status.idle":"2022-08-11T19:22:00.636173Z","shell.execute_reply.started":"2022-08-11T19:22:00.351459Z","shell.execute_reply":"2022-08-11T19:22:00.634585Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Lasso test set score\ntest_scores(lasso_mod, X_test, y_test)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:22:13.406215Z","iopub.execute_input":"2022-08-11T19:22:13.406670Z","iopub.status.idle":"2022-08-11T19:22:13.420467Z","shell.execute_reply.started":"2022-08-11T19:22:13.406632Z","shell.execute_reply":"2022-08-11T19:22:13.417584Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Lasso MSE & RMSE scores\naccuracy_scores(lasso_mod, X_test, y_test)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:22:23.215920Z","iopub.execute_input":"2022-08-11T19:22:23.216306Z","iopub.status.idle":"2022-08-11T19:22:23.233024Z","shell.execute_reply.started":"2022-08-11T19:22:23.216275Z","shell.execute_reply":"2022-08-11T19:22:23.231265Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Collecting Lasso coefficients\ndf_lasso_coef = pd.DataFrame(lasso_mod.coef_, index=X_train.columns,\n                       columns=['Coefficients'])\ndf_lasso_coef['Coef_abs'] = df_lasso_coef.Coefficients.abs()","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:23:58.018193Z","iopub.execute_input":"2022-08-11T19:23:58.018658Z","iopub.status.idle":"2022-08-11T19:23:58.026751Z","shell.execute_reply.started":"2022-08-11T19:23:58.018622Z","shell.execute_reply":"2022-08-11T19:23:58.025733Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Comparing the coefficients from the three linear models","metadata":{}},{"cell_type":"code","source":"coef_0609 = pd.concat([df_ols_coef['Coefficients'], df_ridge_coef['Coefficients'], df_lasso_coef['Coefficients']])\ncoef_0609 = pd.DataFrame(coef_0609)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:24:58.815303Z","iopub.execute_input":"2022-08-11T19:24:58.815908Z","iopub.status.idle":"2022-08-11T19:24:58.824656Z","shell.execute_reply.started":"2022-08-11T19:24:58.815844Z","shell.execute_reply":"2022-08-11T19:24:58.823456Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"coef_0609.reset_index(level=0, inplace=True)\ncoef_0609.columns = ['variable', 'coefficient']\ncoef_0609","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:25:08.486058Z","iopub.execute_input":"2022-08-11T19:25:08.486479Z","iopub.status.idle":"2022-08-11T19:25:08.508872Z","shell.execute_reply.started":"2022-08-11T19:25:08.486449Z","shell.execute_reply":"2022-08-11T19:25:08.507565Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"coef_0609.loc[0:26, \"model\"] = \"ols\"\ncoef_0609.loc[27:53, \"model\"] = \"ridge\"\ncoef_0609.loc[54:80, \"model\"] = \"lasso\"","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:25:32.555968Z","iopub.execute_input":"2022-08-11T19:25:32.556413Z","iopub.status.idle":"2022-08-11T19:25:32.565181Z","shell.execute_reply.started":"2022-08-11T19:25:32.556375Z","shell.execute_reply":"2022-08-11T19:25:32.564081Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"w = sns.catplot(x='variable', y='coefficient', hue='model', data=coef_0609, kind='bar', height=8, aspect=2)\n# set rotation\nw.set_xticklabels(rotation=90)\n\nplt.title('Coefficients of various models on 2006-2009 training set', fontsize=20)\nplt.xlabel(\"Variables\", size=16)\nplt.ylabel(\"Coefficients\", size=16)\nplt.legend(loc=\"upper right\", ncol=2, fontsize=12)\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:25:47.361814Z","iopub.execute_input":"2022-08-11T19:25:47.362215Z","iopub.status.idle":"2022-08-11T19:25:48.658143Z","shell.execute_reply.started":"2022-08-11T19:25:47.362181Z","shell.execute_reply":"2022-08-11T19:25:48.656946Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"The coefficients show that the six more important variables in terms of their impact size on the target are: GrLivArea, OverallQual, Age, OverCond, BaseLivArea, and Location.\n\nBoth the coefficients and R-squared of all three linear models appear very stable (so low variance) across the OLS, Ridge and Lasso models. The R-squared is approximately 0.91-0.92 across all three models, and across the training and test sets too. Moreover, consistent RMSE of approximately 0.1057-0.1059 across the three models.","metadata":{}},{"cell_type":"markdown","source":"### Testing with 2010 Holdout Data","metadata":{}},{"cell_type":"code","source":"y_train = df_0609['LnSalePrice']\ny_lnSP_2010 = df_2010['LnSalePrice']","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:27:09.329375Z","iopub.execute_input":"2022-08-11T19:27:09.329824Z","iopub.status.idle":"2022-08-11T19:27:09.336795Z","shell.execute_reply.started":"2022-08-11T19:27:09.329790Z","shell.execute_reply":"2022-08-11T19:27:09.335037Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"X_2010 = df_2010.drop(['SalePrice', 'LnSalePrice'], axis=1)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:27:20.738509Z","iopub.execute_input":"2022-08-11T19:27:20.739037Z","iopub.status.idle":"2022-08-11T19:27:20.746591Z","shell.execute_reply.started":"2022-08-11T19:27:20.738991Z","shell.execute_reply":"2022-08-11T19:27:20.745211Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Using the full 2006-2009 data to train the model\nX_train = df_0609.drop(['SalePrice', 'LnSalePrice'], axis=1)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:27:28.304745Z","iopub.execute_input":"2022-08-11T19:27:28.305188Z","iopub.status.idle":"2022-08-11T19:27:28.312150Z","shell.execute_reply.started":"2022-08-11T19:27:28.305152Z","shell.execute_reply":"2022-08-11T19:27:28.310955Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Setting 2010 data as the test data\nscaler = StandardScaler()\nX_train = pd.DataFrame(scaler.fit_transform(X_train), columns=X_train.columns)\nX_test = pd.DataFrame(scaler.transform(X_2010), columns=X_2010.columns)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:27:42.534339Z","iopub.execute_input":"2022-08-11T19:27:42.534748Z","iopub.status.idle":"2022-08-11T19:27:42.549948Z","shell.execute_reply.started":"2022-08-11T19:27:42.534709Z","shell.execute_reply":"2022-08-11T19:27:42.548728Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"ols.fit(X_train, y_train) ","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:27:55.472808Z","iopub.execute_input":"2022-08-11T19:27:55.473232Z","iopub.status.idle":"2022-08-11T19:27:55.490962Z","shell.execute_reply.started":"2022-08-11T19:27:55.473198Z","shell.execute_reply":"2022-08-11T19:27:55.489380Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# OLS training set scores, including CV scores\ntrain_scores(ols, X_train, y_train)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:28:05.025719Z","iopub.execute_input":"2022-08-11T19:28:05.026549Z","iopub.status.idle":"2022-08-11T19:28:05.107324Z","shell.execute_reply.started":"2022-08-11T19:28:05.026495Z","shell.execute_reply":"2022-08-11T19:28:05.105933Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# OLS test set score for 2010 holdout set\ntest_scores(ols, X_test, y_lnSP_2010)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:28:19.630556Z","iopub.execute_input":"2022-08-11T19:28:19.631146Z","iopub.status.idle":"2022-08-11T19:28:19.640623Z","shell.execute_reply.started":"2022-08-11T19:28:19.631096Z","shell.execute_reply":"2022-08-11T19:28:19.639690Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# OLS MSE & RMSE scores for 2010 holdout set\naccuracy_scores(ols, X_test, y_lnSP_2010)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:28:30.407074Z","iopub.execute_input":"2022-08-11T19:28:30.407951Z","iopub.status.idle":"2022-08-11T19:28:30.416330Z","shell.execute_reply.started":"2022-08-11T19:28:30.407908Z","shell.execute_reply":"2022-08-11T19:28:30.415118Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"predictions_train = ols.predict(X_train)\npredictions_test = ols.predict(X_test)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:28:43.125323Z","iopub.execute_input":"2022-08-11T19:28:43.125902Z","iopub.status.idle":"2022-08-11T19:28:43.151651Z","shell.execute_reply.started":"2022-08-11T19:28:43.125830Z","shell.execute_reply":"2022-08-11T19:28:43.149134Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Descriptive statistics of model residuals\nols_residuals_2010 = (y_lnSP_2010 - predictions_test)\nols_residuals_2010.describe()","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:28:55.702084Z","iopub.execute_input":"2022-08-11T19:28:55.702492Z","iopub.status.idle":"2022-08-11T19:28:55.715583Z","shell.execute_reply.started":"2022-08-11T19:28:55.702460Z","shell.execute_reply":"2022-08-11T19:28:55.714523Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"print(\"Skew:\", ols_residuals_2010.skew())\nprint(\"Kurtosis:\", ols_residuals_2010.kurtosis())\nstat, p = shapiro(ols_residuals_2010)\nprint('Shapiro-Wilk test on normality=%.3f, p=%.3f' % (stat, p))","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:29:13.834445Z","iopub.execute_input":"2022-08-11T19:29:13.835806Z","iopub.status.idle":"2022-08-11T19:29:13.843523Z","shell.execute_reply.started":"2022-08-11T19:29:13.835755Z","shell.execute_reply":"2022-08-11T19:29:13.842604Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"resid_histogram(ols, X_test, y_lnSP_2010, period='2010 holdout data')","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:29:37.851084Z","iopub.execute_input":"2022-08-11T19:29:37.851536Z","iopub.status.idle":"2022-08-11T19:29:38.312975Z","shell.execute_reply.started":"2022-08-11T19:29:37.851500Z","shell.execute_reply":"2022-08-11T19:29:38.311777Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"fig, ax = plt.subplots(ncols=2, figsize=(15, 6))\nax[0].hist(ols_residuals_2010, density=True, bins=30, color='indianred')\nax[0].set_title('OLS Residuals, 2010 test data', fontsize=14)\nax[1].scatter(ols_residuals_2010, predictions_test, color='midnightblue')\nax[1].set_title('Residuals vs Predicted Target, 2010 test data', fontsize=14)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:29:53.107425Z","iopub.execute_input":"2022-08-11T19:29:53.108488Z","iopub.status.idle":"2022-08-11T19:29:53.709130Z","shell.execute_reply.started":"2022-08-11T19:29:53.108441Z","shell.execute_reply":"2022-08-11T19:29:53.707994Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"from scipy import stats\nstats.probplot(ols_residuals_2010, dist=\"norm\", plot=plt)\nplt.title(\"Quantile-Quantile Plot, 2010 holdout data\")\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:30:07.854839Z","iopub.execute_input":"2022-08-11T19:30:07.855575Z","iopub.status.idle":"2022-08-11T19:30:08.138814Z","shell.execute_reply.started":"2022-08-11T19:30:07.855530Z","shell.execute_reply":"2022-08-11T19:30:08.137536Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Ridge & Lasso regressions on 2010 data","metadata":{}},{"cell_type":"code","source":"ridge_mod = RidgeCV(alphas=np.logspace(-4, 4, 10), cv=5)\nridge_mod.fit(X_train, y_train)\n\nprint('Best Ridge alpha:', ridge_mod.alpha_)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:30:35.204323Z","iopub.execute_input":"2022-08-11T19:30:35.204736Z","iopub.status.idle":"2022-08-11T19:30:35.906010Z","shell.execute_reply.started":"2022-08-11T19:30:35.204701Z","shell.execute_reply":"2022-08-11T19:30:35.904019Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"train_scores(ridge_mod, X_train, y_train)\ntest_scores(ridge_mod, X_test, y_lnSP_2010)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:30:46.525999Z","iopub.execute_input":"2022-08-11T19:30:46.527186Z","iopub.status.idle":"2022-08-11T19:30:49.270232Z","shell.execute_reply.started":"2022-08-11T19:30:46.527142Z","shell.execute_reply":"2022-08-11T19:30:49.268552Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"lasso_mod = LassoCV(alphas=np.logspace(-4, 4, 10), cv=5)\nlasso_mod.fit(X_train, y_train)\n\nprint('Best Lasso alpha:', lasso_mod.alpha_)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:31:02.263271Z","iopub.execute_input":"2022-08-11T19:31:02.263728Z","iopub.status.idle":"2022-08-11T19:31:02.311621Z","shell.execute_reply.started":"2022-08-11T19:31:02.263694Z","shell.execute_reply":"2022-08-11T19:31:02.310031Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"train_scores(lasso_mod, X_train, y_train)\ntest_scores(lasso_mod, X_test, y_lnSP_2010)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:31:12.267176Z","iopub.execute_input":"2022-08-11T19:31:12.267637Z","iopub.status.idle":"2022-08-11T19:31:12.571924Z","shell.execute_reply.started":"2022-08-11T19:31:12.267599Z","shell.execute_reply":"2022-08-11T19:31:12.570143Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Regression model using full 2006-2010 data","metadata":{}},{"cell_type":"code","source":"y_lnSP = df['LnSalePrice']","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:32:25.282338Z","iopub.execute_input":"2022-08-11T19:32:25.282771Z","iopub.status.idle":"2022-08-11T19:32:25.288825Z","shell.execute_reply.started":"2022-08-11T19:32:25.282735Z","shell.execute_reply":"2022-08-11T19:32:25.287690Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"X_fin = df.drop(['SalePrice', 'LnSalePrice'], axis=1)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:32:32.461099Z","iopub.execute_input":"2022-08-11T19:32:32.461825Z","iopub.status.idle":"2022-08-11T19:32:32.468361Z","shell.execute_reply.started":"2022-08-11T19:32:32.461788Z","shell.execute_reply":"2022-08-11T19:32:32.467388Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"X_fin.shape","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:32:39.094565Z","iopub.execute_input":"2022-08-11T19:32:39.095086Z","iopub.status.idle":"2022-08-11T19:32:39.103541Z","shell.execute_reply.started":"2022-08-11T19:32:39.095038Z","shell.execute_reply":"2022-08-11T19:32:39.102621Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"X_fin = pd.DataFrame(scaler.fit_transform(X_fin), columns=X_fin.columns)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:32:55.213917Z","iopub.execute_input":"2022-08-11T19:32:55.214333Z","iopub.status.idle":"2022-08-11T19:32:55.225465Z","shell.execute_reply.started":"2022-08-11T19:32:55.214300Z","shell.execute_reply":"2022-08-11T19:32:55.224242Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"ols.fit(X_fin, y_lnSP)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:33:07.043648Z","iopub.execute_input":"2022-08-11T19:33:07.044110Z","iopub.status.idle":"2022-08-11T19:33:07.062494Z","shell.execute_reply.started":"2022-08-11T19:33:07.044072Z","shell.execute_reply":"2022-08-11T19:33:07.060680Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Model scores on the full 2006-2010 data\ntest_scores(ols, X_fin, y_lnSP)\naccuracy_scores(ols, X_fin, y_lnSP)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:33:16.835730Z","iopub.execute_input":"2022-08-11T19:33:16.836213Z","iopub.status.idle":"2022-08-11T19:33:16.860011Z","shell.execute_reply.started":"2022-08-11T19:33:16.836176Z","shell.execute_reply":"2022-08-11T19:33:16.858449Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_LnSP_coef = pd.DataFrame(ols.coef_, index=X_fin.columns,\n                       columns=['Coefficients'])\ndf_LnSP_coef['Coef_abs'] = df_LnSP_coef.Coefficients.abs()","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:33:27.545102Z","iopub.execute_input":"2022-08-11T19:33:27.546172Z","iopub.status.idle":"2022-08-11T19:33:27.553664Z","shell.execute_reply.started":"2022-08-11T19:33:27.546129Z","shell.execute_reply":"2022-08-11T19:33:27.552169Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Descriptive statistics of the residuals\npredictions = ols.predict(X_fin)\nerror_term = (y_lnSP - predictions)\nerror_term.describe()","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:33:37.914467Z","iopub.execute_input":"2022-08-11T19:33:37.915083Z","iopub.status.idle":"2022-08-11T19:33:37.957641Z","shell.execute_reply.started":"2022-08-11T19:33:37.915040Z","shell.execute_reply":"2022-08-11T19:33:37.955931Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Residuals are approximately normal, given the skew and kurtosis\nprint(\"Skew:\", error_term.skew())\nprint(\"Kurtosis:\", error_term.kurtosis())\nstat, p = shapiro(error_term)\nprint('Shapiro-Wilk test on normality=%.3f, p=%.3f' % (stat, p))","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:33:49.031092Z","iopub.execute_input":"2022-08-11T19:33:49.031513Z","iopub.status.idle":"2022-08-11T19:33:49.040113Z","shell.execute_reply.started":"2022-08-11T19:33:49.031477Z","shell.execute_reply":"2022-08-11T19:33:49.038913Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"resid_histogram(ols, X_fin, y_lnSP, period='full 2006-2010 data')","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:34:00.644602Z","iopub.execute_input":"2022-08-11T19:34:00.645530Z","iopub.status.idle":"2022-08-11T19:34:01.158174Z","shell.execute_reply.started":"2022-08-11T19:34:00.645492Z","shell.execute_reply":"2022-08-11T19:34:01.156810Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"from scipy import stats\nstats.probplot(error_term, dist=\"norm\", plot=plt)\nplt.title(\"Quantile-Quantile Plot, 2006-2010 data\")\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:34:15.134340Z","iopub.execute_input":"2022-08-11T19:34:15.134735Z","iopub.status.idle":"2022-08-11T19:34:15.424891Z","shell.execute_reply.started":"2022-08-11T19:34:15.134704Z","shell.execute_reply":"2022-08-11T19:34:15.423690Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Statistical inference and hypothesis testing","metadata":{}},{"cell_type":"code","source":"# Using the full data matrix to get the unstandardised coefficients\nX_fin = df.drop(['SalePrice', 'LnSalePrice'], axis=1)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:34:37.501620Z","iopub.execute_input":"2022-08-11T19:34:37.502067Z","iopub.status.idle":"2022-08-11T19:34:37.508411Z","shell.execute_reply.started":"2022-08-11T19:34:37.502029Z","shell.execute_reply":"2022-08-11T19:34:37.507350Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Need to add a constant for the Statsmodel OLS model\nX_sm = sm.add_constant(X_fin)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:34:46.413753Z","iopub.execute_input":"2022-08-11T19:34:46.414165Z","iopub.status.idle":"2022-08-11T19:34:46.431679Z","shell.execute_reply.started":"2022-08-11T19:34:46.414131Z","shell.execute_reply":"2022-08-11T19:34:46.430814Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Model using full 2006-2010 data with unstandardised values, featuring individual coefficients with p-values\nmodel = sm.OLS(y_lnSP, X_sm)\nresults = model.fit()\nprint(results.summary())","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:34:57.014957Z","iopub.execute_input":"2022-08-11T19:34:57.015494Z","iopub.status.idle":"2022-08-11T19:34:57.070944Z","shell.execute_reply.started":"2022-08-11T19:34:57.015445Z","shell.execute_reply":"2022-08-11T19:34:57.069290Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"The **majority of the variables** in the OLS model are **statistically significant** at the 5% level, as indicated by the individual t-statistics and p-values. The exceptions are 'RoadRail' and 'LowQualFinSF'.","metadata":{}},{"cell_type":"code","source":"ols.fit(X_fin, y_lnSP)\nunscaled_coef = ols.coef_","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:35:48.339951Z","iopub.execute_input":"2022-08-11T19:35:48.340425Z","iopub.status.idle":"2022-08-11T19:35:48.358689Z","shell.execute_reply.started":"2022-08-11T19:35:48.340384Z","shell.execute_reply":"2022-08-11T19:35:48.356522Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import math\n\ntransformed_coef = []\nfor i in unscaled_coef:\n    j = math.exp(i)\n    transformed_coef.append(j)\nprint(transformed_coef)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:35:57.376461Z","iopub.execute_input":"2022-08-11T19:35:57.377017Z","iopub.status.idle":"2022-08-11T19:35:57.385036Z","shell.execute_reply.started":"2022-08-11T19:35:57.376969Z","shell.execute_reply":"2022-08-11T19:35:57.383740Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"coef_effect = [(i - 1)*df.SalePrice.mean() for i in transformed_coef]","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:36:09.293932Z","iopub.execute_input":"2022-08-11T19:36:09.294391Z","iopub.status.idle":"2022-08-11T19:36:09.304402Z","shell.execute_reply.started":"2022-08-11T19:36:09.294348Z","shell.execute_reply":"2022-08-11T19:36:09.303402Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"var_impact = pd.DataFrame(data=[X_fin.columns, coef_effect]).T","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:36:17.173813Z","iopub.execute_input":"2022-08-11T19:36:17.174309Z","iopub.status.idle":"2022-08-11T19:36:17.184015Z","shell.execute_reply.started":"2022-08-11T19:36:17.174268Z","shell.execute_reply":"2022-08-11T19:36:17.182690Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"var_impact.columns = [\"variable\", \"1-unit change\"]","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:36:27.626505Z","iopub.execute_input":"2022-08-11T19:36:27.626975Z","iopub.status.idle":"2022-08-11T19:36:27.632879Z","shell.execute_reply.started":"2022-08-11T19:36:27.626936Z","shell.execute_reply":"2022-08-11T19:36:27.631962Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"var_impact","metadata":{"execution":{"iopub.status.busy":"2022-08-11T19:36:34.912846Z","iopub.execute_input":"2022-08-11T19:36:34.913356Z","iopub.status.idle":"2022-08-11T19:36:34.925787Z","shell.execute_reply.started":"2022-08-11T19:36:34.913316Z","shell.execute_reply":"2022-08-11T19:36:34.924635Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### ","metadata":{}}]}