<jupyter_start><jupyter_text>real-estate  AbstractThis report goes explores a highly dimensional dataset on properties to find a regression model for the sale price. During the first analysis stage, the data dimensions are reduced from 82 to ~30. Further reductions are achieved by closer analysis.It was found that methods such as PCA (on the continuous part of the data), MCA (on the categorical parts of the data) and factor analysis did not result in explainable attributes, and thus could not be used to reduce the dimensions much further.In the end, a Ridge regression model was chosen, with $\lambda = 8$. This model provided a fairly low MAE = 0.15, but also a fairly low R^2 $\approx$ 0.48. The highlight of this project has been in it's use of the pandas_profiling library to quickly remove attributes that are deemed uninportant.  Scope and ObjectivesThe purpose of this notebook is to develop a regression model for the sale price of houses. Note that the author, at this point in time, is otherwise occupied with other work and has thus been unable to perform a comprehensive analysis for this project.Therefore, to accelerate the process, the author placed a great amount of importance into using an external library, pandas_profiling, to accelerate the process of reducing dimensions.Justifications for removing all of the variables will not be provided, but some specific examples will be.The project objectives are to:1) develop a regression model that could predict the sale price of property2) reduce the dimensions of the data3) explore specifically the concept of PCA, MCA and factor analysis for dimensionality reduction  Libraries<jupyter_code># ICDSS
from PythonFiles import reader as r

# mathematical
import numpy as np

# data storage / analysis
import pandas as pd
from pandas_profiling import ProfileReport

# dimensional reduction / statistical analysis
import prince
from sklearn.decomposition import PCA
from scipy.stats import skew as skewness

# plotting
import matplotlib.pyplot as plt
import seaborn as sns
from scipy.stats import probplot

# preprocessing
from sklearn.preprocessing import StandardScaler
from sklearn.model_selection import train_test_split

# modelling
from sklearn.linear_model import LinearRegression
from sklearn.linear_model import Ridge

# evaluation
from sklearn.metrics import r2_score
from sklearn.metrics import mean_absolute_error as MAE<jupyter_output><empty_output><jupyter_text>Exploration - Reducing Dimensions<jupyter_code>path = 'Data/data.csv'
df = r.get_housing_data(path) # unnecessary change to your function, just did this out of curiosity
df.head()<jupyter_output><empty_output><jupyter_text>**NOTE:** running the next 3 cells (already commneted) will take a considerable amount of time. The outcome of these 3 cells can be found in 'Reporting/profile.html'.<jupyter_code>#profile = df.profile_report()
#profile

#profile.to_file('Reporting/profile.html')
df_temp = df.copy()

# drop attributes which contain > 30% missing values
df_temp.drop(['PoolQC','Fence','Alley','FireplaceQu','MiscFeature'], axis = 1, inplace = True)

# drop attributes that are redundant / meaningless based on our knowledge
df_temp.drop(['Order','PID','SaleType'], axis = 1, inplace = True)

# drop attributes that have very little correlation with everything (based on the correlation graph generated)
df_temp.drop(['3SsnPorch', 'ScreenPorch','PoolArea',
              'MiscVal','MoSold','YrSold','MSZoning',
              'MSSubClass','LotFrontage','LotArea',
              'LotShape','LotConfig','BldgType',
              'HouseStyle'], axis = 1, inplace = True)

# drop attributes that are highly skewed with weak correlations
# Examples: 99.5% of 'street' is 'Pave'

df_temp.drop(['Street','LandContour','Utilities','LandSlope',
              'Condition1','Condition2','RoofMatl','BsmtCond','PavedDrive',
              'KitchenAbvGr','Functional','EnclosedPorch','CentralAir',
              'GarageCond','GarageQual','BsmtHalfBath','LowQualFinSF','Electrical',
              'LandContour','LandSlope','ExterCond'], axis = 1, inplace = True)

# drop values that may mislead our analysis (this includes things that are incorrectly
# stored, for instance: year remod/add = year constructed if no remodelling as been done)
df_temp.drop(['YearRemod/Add','GarageCars'],axis = 1, inplace = True)

df_temp.head()
# replace missing values with median for numerical, mode for categorical
# note that for these columns, the percentage of missing values is < 6% (except for LotFrontage at 16.7%)
for column in df_temp.columns:
    if df_temp[column].dtypes in ['int64', 'float64']:
        df_temp[column] = df_temp[column].fillna(df_temp[column].median())
    else:
        df_temp[column] = df_temp[column].fillna(df_temp[column].mode()[0])
# verify that data types are correct (some nominal values may be read as 'int')
df_temp.info()
# change variables that have incorrect data types
# example: the year a building was built should not be continuous data
df_temp[['GarageYrBlt','YearBuilt']] = df_temp[['GarageYrBlt','YearBuilt']].astype(str)
#profile_temp = df_temp.profile_report()

#profile_temp
#profile_temp.to_file('Reporting/profile2.html')
# further cancellations
df_temp2 = df_temp.copy()
df_temp2.drop(['BsmtFinType2','BsmtFinSF2','Heating',
              'BsmtFullBath','HalfBath','BedroomAbvGr',
               'WoodDeckSF','OpenPorchSF','RoofStyle'], axis =1, inplace = True)

df_temp2.head()<jupyter_output><empty_output><jupyter_text>For the time being, reducing the initial 82 attributes to 30 has been good. Since the remaining values are somewhat intercorrelated (highly, at least in comparison to the previous varriables), it might be wise to perform some other dimensionality reduction methods, such as PCA.<jupyter_code>(df_temp2.dtypes == 'object').value_counts(True)<jupyter_output><empty_output><jupyter_text>Given that almost half of the attributes are actually categorical variables, this makes it difficult to justify the use of PCA analysis (one-hot encoded variables will not fix the problem, due to the way PCA operates). As a result, a factor analysis (this *sort of* combines MCA and PCA) will be used to reduce the dimensions even further.  Exploration - Reducing dimensions<jupyter_code># normalize non-object data
scaler = StandardScaler()

df_temp3 = df_temp2.copy()
df_temp3[df_temp3.columns[(df_temp3.dtypes != 'object')]] =\
scaler.fit_transform(df_temp3[df_temp3.columns[(df_temp3.dtypes != 'object')]])

df_temp3.head()
famd = prince.FAMD(
    n_components = 25,
    n_iter = 10
)
famd = famd.fit(df_temp3.drop('SalePrice', axis='columns'))

plt.plot([i for i in range(1,26)],famd.eigenvalues_)<jupyter_output><empty_output><jupyter_text>Based on the above, it suffices to say that Factor Analysis does not in fact provide suitable results (this is because the eigenvectors are all below 1). Let's try to see if we can reduce either the dimensions of the categorical data using MCA, or that of the continuous data using PCA<jupyter_code>pca = PCA(
    n_components = 13
)

df_numerical = df_temp3[df_temp3.columns[df_temp3.dtypes != 'object']]
df_categorical = df_temp3[df_temp3.columns[df_temp3.dtypes == 'object']]

Y= pca.fit_transform(np.asarray(df_numerical.drop('SalePrice', axis = 1)))
pca.explained_variance_ratio_.cumsum()
# choosing n_components = 8 will retain 94% of the variance, while reducing our
# dimensions by 6
pca = PCA(
    n_components = 8
)

df_numerical = pca.fit_transform(np.asarray(df_numerical.drop('SalePrice', axis =1)))

fig = plt.figure(figsize = (16,16))
for i in range(8):
    fig.add_subplot(4,2,i+1)
    plt.plot(df_numerical[:,i], df_temp3['SalePrice'],'.')

# not exactly convinced by the outcome of the PCA<jupyter_output><empty_output><jupyter_text>This method of reducing the dimensions of the variables is not convincing, and in fact is problematic: we've essentially *lost* the meaning of the variables. So, even if the PCA (or factor analysis for that matter) gave us 4-5 variables in total, the models produced would not be insightful since we'd struggle to interpret them.  Modelling<jupyter_code># convert categorical varialbes to dummy ones
df_numerical = df_temp3[df_temp3.columns[df_temp3.dtypes != 'object']]
df_dummies = pd.get_dummies(data = df_categorical)

df_temp4 = pd.concat([df_numerical, df_dummies], axis = 1)
df_temp4.shape<jupyter_output><empty_output><jupyter_text>Since we are interested in modelling the sale price, we are not exactly interested in understanding why outliers occur. As such, it might be worth removing them. To do this, we need to examine how skewed out datasets are. If they are extremely skewed, then we may benefit from performing a de-skewing transformation before removing outliers (as this would make the outlier removal symmetric).<jupyter_code>skewnesses = []
for i in df_numerical.columns:
    skew = skewness(df_numerical[i])
    if skew > 1:
        print(i,skew)
        df_numerical[i+'_sqrt'] = np.sqrt(np.abs(df_numerical[i]))
        df_numerical.drop(i,axis =1, inplace = True)
        skew = skewness(df_numerical[i+'_sqrt'])
        print(skew)
    skewnesses.append(skew)

print(skewnesses)
# remove outliers

fig, ax = plt.subplots(figsize = (16,10))
ax = sns.boxplot(ax = ax,data=df_numerical, orient="v", color="skyblue")
ax.set_xticklabels(ax.get_xticklabels(), rotation=20)
ax<jupyter_output><empty_output><jupyter_text>Seeing as the MasVnrArea was mostly comprised of outliers, it was decided to explore this a bit further.It's worth noting that this attribute still remained higlhy skewed, even after transformation. It is mostly comprised of zeroes, however it seems to be correlated with house SalePrice.<jupyter_code># examine the box plot for this attribute before transfrmation
fig, ax = plt.subplots(figsize = (16,10))
ax = sns.boxplot(ax = ax,x=df_temp4['MasVnrArea'], orient="h", color="skyblue")
ax
fig, ax = plt.subplots(figsize = (16,10))
ax = sns.boxplot(ax = ax,x=np.sqrt(np.abs(df_temp4['MasVnrArea'])), orient="h", color="skyblue")
ax<jupyter_output><empty_output><jupyter_text>Upon further inspection, it makes sense that this is, somewhat expected. That said, it would be illogical to consdier anything over UQ + 1.25IQR and anything lower than LQ - 1.25IQR as outliers, since thst would essentially kill the dataset.An alternative approach would be to define more stringent conditions for this particular attribute, so for instance UQ + $\gamma$IQR, where $\gamma$ is set to a high value. However, even this is arbitrary, and we have no knowledge of what value to set. With that in mind, it was decided to dock this attribute off.<jupyter_code>df_numerical.drop('MasVnrArea_sqrt', axis =1 , inplace= True)
# remove  values that exceed UQ + 1.25IQR and those below LQ - 1.25IQR

for i in df_numerical.columns:
    LQ = df_numerical[i].quantile(0.25)
    UQ = df_numerical[i].quantile(0.75)
    IQR = UQ - LQ
    df_numerical = df_numerical[(df_numerical[i] < 1.25*IQR + UQ) &
                                (df_numerical[i] > LQ - 1.25*IQR)
                               ]
df_numerical.shape[0]/2930
# dataset has been decreased by 30%

# now only store the indices that are relevant

df_temp5 = df_numerical.copy()

df_temp5 = pd.concat([df_temp5, df_dummies.iloc[df_temp5.index]], axis = 1)

df_temp5.to_csv('Data/data_reduced.csv')

df_temp5.shape # correctly concatenated
X = df_temp5[df_temp5.columns[df_temp5.columns != 'SalePrice_sqrt']]
y = df_temp5['SalePrice_sqrt']
X_train, X_test, y_train, y_test = train_test_split(
    X, y, test_size=0.3)

model = LinearRegression()
model.fit(X_train,y_train)
y_hat = model.predict(X_test)

sns.regplot(y_test, y_hat)<jupyter_output>/Users/yousefnami/python_environments/ICDSS_Challenges/real-estate/lib/python3.7/site-packages/seaborn/_decorators.py:43: FutureWarning: Pass the following variables as keyword args: x, y. From version 0.12, the only valid positional argument will be `data`, and passing other arguments without an explicit keyword will result in an error or misinterpretation.
  FutureWarning<jupyter_text>Clearly, using a linear regression does not seem appropriate. We see that the model is punished by very large errors. This suggests that our method of removing outliers (by removing values exceeding UQ + 1.25IQR and those lower than LQ - 1.25IQR) has not been entirely successful. This calls for using more sophisticated methods of regression.<jupyter_code>model = Ridge()
model.fit(X_train,y_train)
y_hat = model.predict(X_test)
sns.regplot(y_test, y_hat)
plt.show()<jupyter_output>/Users/yousefnami/python_environments/ICDSS_Challenges/real-estate/lib/python3.7/site-packages/seaborn/_decorators.py:43: FutureWarning: Pass the following variables as keyword args: x, y. From version 0.12, the only valid positional argument will be `data`, and passing other arguments without an explicit keyword will result in an error or misinterpretation.
  FutureWarning<jupyter_text>Note: here we invoked the use of the Ridge regression technique, because of it's regularisation parameter ($\lambda$) that decreases the effect of high errors. It's important to note that changing this value will affect the bias and variance of the model. Increasing it will lead to a high balance, but low variance, whereas a decrease causes the opposite. Let's plot to see what the most suitable value is.<jupyter_code>r2 = []
for i in range(1, 20):
    model = Ridge(alpha = i)
    model.fit(X_train,y_train)
    y_hat = model.predict(X_test)
    r2.append(r2_score(y_test,y_hat))
    
    
l1 = plt.plot([i for i in range(1,20)], r2)
l2 = plt.plot(r2.index(max(r2)),max(r2),'rx')
plt.legend(l2,['max $R^2$ at $\\alpha = 8$'])<jupyter_output><empty_output><jupyter_text>Evaluating model performance As seen from above, the $R^2$ of the model is fairly low. However, this could be due to high variance in the data. We can get a better understanding of how good the model is by exploring some other error metrics.<jupyter_code>model = Ridge(alpha = 8)
model.fit(X_train,y_train)
y_hat = model.predict(X_test)
MAE(y_test, y_hat)<jupyter_output><empty_output><jupyter_text>It seems upon closer inspection that the model has not in fact performed terribly. The low R^2 value is merely suggestive of the fact that the data has high variance.<jupyter_code># check the data for heteroscedasticity
error = y_test - y_hat
plt.plot(error,'.')
# check errors for normality
fig, ax = plt.subplots(figsize=(16,8))
probplot(error, plot = ax)
plt.show()<jupyter_output><empty_output><jupyter_text>From the two plots above, it is clear that the assumption of normal errors is valid.<jupyter_code># explore the ridge coefficients
coefficients = pd.DataFrame(data = {'Attribute':X.columns, 'Coefficient':model.coef_})
plt.hist(coefficients['Coefficient'])<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.09329917756505325, 'fraction_numerical': 0.014089254415531886, 'mean_word_length': 5.201923076923077, 'pattern_counts': {'":': 0, '<': 61, '<?xml version=': 0, '>': 62, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 19, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '15', 'dir': 'jupyter-structured-clean-dedup', 'id': '553974', 'license': 'MIT', 'n_tokens_mistral': 4742, 'n_tokens_neox': 4464, 'n_words': 1585, 'path': 'real-estate/real-estate.ipynb', 'repo_name': 'namiyousef/ICDSS_Challenges2020', 'size': '14834'}