<jupyter_start><jupyter_text>House Price Prediction  Importing Libraries<jupyter_code>import os
import numpy as np
import matplotlib.pyplot as plt
import tensorflow as tf
import pathlib
import pandas as pd
import urllib
import tarfile
import sklearn
import joblib

from scipy.stats import *
from sklearn import *
from sklearn.model_selection import *

print(tf.__version__)
print(sklearn.__version__)<jupyter_output>2.3.0
0.23.2<jupyter_text>Data Acquisition and Visualization<jupyter_code>dataset_url = "https://github.com/ageron/handson-ml2/raw/master/datasets/housing/housing.tgz"
data_dir = keras.utils.get_file('housing', origin = dataset_url, untar = True, extract=True)
data_dir = pathlib.Path(data_dir)
# print(data_dir)
df = pd.read_csv(os.path.join(data_dir.parents[0], 'housing.csv'))
df.head()
df.info()
df['ocean_proximity'].value_counts()
df.hist(bins = 50, figsize = (15,12))
plt.show()<jupyter_output><empty_output><jupyter_text>Data Preprocessing<jupyter_code>df['income_cat'] = pd.cut(df['median_income'], bins = [0, 1.5, 3, 4.5, 6, np.inf], labels = [1,2,3,4,5])
df['income_cat'].hist()
plt.show()<jupyter_output><empty_output><jupyter_text>**Stratified sampling** divides the entire dataset into homogeneous subgroups. This reduces the sampling bias and ensures the test set perserves the statistical distribution of the data (representative).Shouldn not make too many strata (sufficient samples in each strata), otherwise the importance estimated most likely will be biased.<jupyter_code>for train_idx, test_idx in StratifiedShuffleSplit(n_splits=1, test_size=0.2, random_state=42).split(df, df['income_cat']):
    train_ds = df.loc[train_idx]
    test_ds = df.loc[test_idx]
def proportion(data):
    return data['income_cat'].value_counts()/len(data)

train_ds_orig, test_ds_orig = train_test_split(df, test_size = 0.2, random_state = 42)

compare_df = pd.DataFrame({
    'Overall': proportion(df),
    'Stratified': proportion(train_ds),
    'Random': proportion(test_ds_orig)
}).sort_index()

compare_df['Rand Error %'] = compare_df['Random'] / compare_df['Overall'] - 1
compare_df['Stra Error %'] = compare_df['Stratified'] / compare_df['Overall'] - 1
compare_df
for set_ in (train_ds, test_ds):
    set_.drop('income_cat', axis = 1, inplace = True)<jupyter_output><empty_output><jupyter_text>Data Investigation and Analysis<jupyter_code>data = train_ds.copy()
data.plot(kind = 'scatter', x = 'longitude', y = 'latitude', alpha = 0.5,
         s = data['population']/100, label = 'Population',
         c = data['median_house_value'], cmap = plt.get_cmap("jet"), figsize = (9,6), colorbar = True)
plt.legend()
plt.show()
corr_matrix = data.corr()
print(corr_matrix['median_house_value'].sort_values(ascending = False))
attributes = ['median_house_value', 'median_income', 'total_rooms', 'housing_median_age']
pd.plotting.scatter_matrix(data[attributes], figsize = (10,8))
data.plot(kind = 'scatter', x = 'median_income', y = 'median_house_value', alpha = 0.2)
plt.show()
data["rooms_per_household"] = data["total_rooms"] / data["households"]
data["bedrooms_per_room"] = data["total_bedrooms"] / data["total_rooms"]
data["population_per_household"] = data["population"] / data["households"]

corr_matrix = data.corr()
corr_matrix['median_house_value'].sort_values(ascending=False)<jupyter_output><empty_output><jupyter_text>It can be noticed that the new features are indeed helpful. For instance, the "rooms_per_household" has a positive correlation (ranks the 3rd place), and it makes sense.  Data Cleaning<jupyter_code># Separate features and labels
train_ds_x = train_ds.drop("median_house_value", axis = 1)
train_ds_y = train_ds["median_house_value"].copy()
# imputer. Using to fill all the N/A cells
train_ds_x_numeric = train_ds_x.drop('ocean_proximity', axis = 1) # Numerical columns

# Do it individually
imputer = sklearn.impute.SimpleImputer(strategy = 'median')
imputer.fit(train_ds_x_numeric)
X = imputer.transform(train_ds_x_numeric)
imputer.statistics_
# Categorical input processing using one-hot-encoder
train_cat = train_ds_x[['ocean_proximity']] # Categorical column

# Do it individually
one_hot_encoder = sklearn.preprocessing.OneHotEncoder()
train_cat_1hot = one_hot_encoder.fit_transform(train_cat)
one_hot_encoder.categories_
# A custom transformer for creating new useful features by combining different features

rooms_ix, bedrooms_ix, population_ix, households_ix = 3, 4, 5, 6 # from the column indices

class CombineAttribsAdder(sklearn.base.BaseEstimator, sklearn.base.TransformerMixin):
    def __init__(self, add_bedrooms_per_room = True): 
        self.add_bedrooms_per_room = add_bedrooms_per_room
    
    def fit(self, X, y=None):
        return self 

    def transform(self, X):
        rooms_per_household = X[:, rooms_ix] / X[:, households_ix]
        population_per_household = X[:, population_ix] / X[:, households_ix]
        if self.add_bedrooms_per_room:
            bedrooms_per_room = X[:, bedrooms_ix] / X[:, rooms_ix]
            return np.c_[X, rooms_per_household, population_per_household, bedrooms_per_room]
        else:
            return np.c_[X, rooms_per_household, population_per_household]
        
attr_adder = CombineAttribsAdder(add_bedrooms_per_room = False)
extra_attribs = attr_adder.transform(train_ds.values)
# Numerical input pipeline
numeric_pipeline = sklearn.pipeline.Pipeline([
    ('imputer', sklearn.impute.SimpleImputer(strategy = 'median')),
    ('attribs_adder', CombineAttribsAdder()),
    ('std_scaler', sklearn.preprocessing.StandardScaler())
])

train_ds_x_num = numeric_pipeline.fit(train_ds_x_numeric)<jupyter_output><empty_output><jupyter_text>All in one using the sklearn.ColummTransformer<jupyter_code>num_attribs = list(train_ds_x_numeric)
cat_attribs = ["ocean_proximity"]

full_pipeline = sklearn.compose.ColumnTransformer([
    ('num', numeric_pipeline, num_attribs),
    ('cat', sklearn.preprocessing.OneHotEncoder(), cat_attribs)
])

train_ds_x = full_pipeline.fit_transform(train_ds)<jupyter_output><empty_output><jupyter_text>Model Building and Selection<jupyter_code>data_pt_input = train_ds_x[:5]
data_pt_label = train_ds_y[:5]<jupyter_output><empty_output><jupyter_text>Linear regression<jupyter_code>lin_reg = linear_model.LinearRegression()
lin_reg.fit(train_ds_x, train_ds_y)
predictions = lin_reg.predict(data_pt_input)
print("predictions\tTrue Labels")
for i, y_hat in enumerate(predictions):
    print('{:.3f}\t{}'.format(y_hat, data_pt_label.to_numpy()[i]))
lin_rmse = np.sqrt(sklearn.metrics.mean_squared_error(data_pt_label, predictions))
print('RMSE: {:.3f}'.format(lin_rmse))<jupyter_output>predictions	True Labels
210644.605	286600.0
317768.807	340600.0
210956.433	196900.0
59218.989	46300.0
189747.558	254500.0
RMSE: 46578.638<jupyter_text>Random Forest<jupyter_code>forest_reg = ensemble.RandomForestRegressor()
forest_reg.fit(train_ds_x, train_ds_y)
predictions = forest_reg.predict(data_pt_input)
print("predictions\tTrue Labels")
for i, y_hat in enumerate(predictions):
    print('{:.3f}\t{}'.format(y_hat, data_pt_label.to_numpy()[i]))
forest_rmse = np.sqrt(metrics.mean_squared_error(data_pt_label, predictions))
print('RMSE: {:.3f}'.format(forest_rmse))<jupyter_output>predictions	True Labels
259186.000	286600.0
328854.000	340600.0
221934.000	196900.0
49396.000	46300.0
238615.000	254500.0
RMSE: 18857.972<jupyter_text>SVM<jupyter_code>svm_reg = svm.SVR(kernel="linear")
svm_reg.fit(train_ds_x, train_ds_y)
predictions = svm_reg.predict(data_pt_input)
print("predictions\tTrue Labels")
for i, y_hat in enumerate(predictions):
    print('{:.3f}\t{}'.format(y_hat, data_pt_label.to_numpy()[i]))
svm_rmse = np.sqrt(metrics.mean_squared_error(data_pt_label, predictions))
print('RMSE: {:.3f}'.format(svm_rmse))<jupyter_output>predictions	True Labels
175421.011	286600.0
190268.746	340600.0
175653.267	196900.0
165556.979	46300.0
183911.428	254500.0
RMSE: 104514.737<jupyter_text>Cross-validation<jupyter_code>def display_score(scores):
    print("Scores: ", scores)
    print('Mean: ', scores.mean())
    print('STD: ', scores.std())
    print("------------------------------------------------------------")
lin_rmse_score = np.sqrt(-cross_val_score(lin_reg, train_ds_x, train_ds_y, scoring="neg_mean_squared_error", cv=10))
forest_rmse_score = np.sqrt(-cross_val_score(forest_reg, train_ds_x, train_ds_y, scoring="neg_mean_squared_error", cv=10))
svm_rmse_score = np.sqrt(-cross_val_score(svm_reg, train_ds_x, train_ds_y, scoring="neg_mean_squared_error", cv=10))
display_score(lin_rmse_score)
display_score(forest_rmse_score)
display_score(svm_rmse_score)<jupyter_output>Scores:  [66782.73843989 66960.118071   70347.95244419 74739.57052552
 68031.13388938 71193.84183426 64969.63056405 68281.61137997
 71552.91566558 67665.10082067]
Mean:  69052.46136345083
STD:  2731.6740017983434
------------------------------------------------------------
Scores:  [49500.78174186 47519.63444107 49866.97178105 52098.64181991
 49805.3641604  53788.9264113  48539.22385948 48193.5496392
 53208.6969429  50258.33818989]
Mean:  50278.012898705565
STD:  2005.092652293203
------------------------------------------------------------
Scores:  [105342.09141998 112489.24624123 110092.35042753 113403.22892482
 110638.90119657 115675.8320024  110703.56887243 114476.89008206
 113756.17971227 111520.1120808 ]
Mean:  111809.84009600841
STD:  2762.393664321567
------------------------------------------------------------<jupyter_text>Fine-tuning When you have no idea what value a hyperparameter should have, a simple approach is to try out consecutive powers of 10, or asmaller number if you want a more fine-grained search  Grid Search<jupyter_code>para_grid = [
    {'n_estimators': [3, 10, 30], 'max_features': [2, 4, 6, 8]},
    {'bootstrap': [False], 'n_estimators': [3, 10], 'max_features': [2, 3, 4]},
]

grid_search = GridSearchCV(forest_reg, para_grid, cv = 5, scoring = 'neg_mean_squared_error', return_train_score = True)
grid_search.fit(train_ds_x, train_ds_y)
grid_search.best_params_
grid_search.best_estimator_
cvres = grid_search.cv_results_
for mean_score, params in zip(cvres["mean_test_score"], cvres["params"]):
    print(np.sqrt(-mean_score), params)
feature_importance = grid_search.best_estimator_.feature_importances_
extra_attribs = ["rooms_per_hhold", "pop_per_hhold", "bedrooms_per_room"]
cat_encoder = full_pipeline.named_transformers_["cat"]
cat_one_hot_attribs = list(cat_encoder.categories_[0])
attributes = num_attribs + extra_attribs + cat_one_hot_attribs
sorted(zip(feature_importance, attributes), reverse=True)
def indices_of_top_k(arr, k):
    return np.sort(np.argpartition(np.array(arr), -k)[-k:])

class TopFeatureSelector(sklearn.base.BaseEstimator, sklearn.base.TransformerMixin):
    def __init__(self, feature_importances, k):
        self.feature_importances = feature_importances
        self.k = k
    def fit(self, X, y=None):
        self.feature_indices_ = indices_of_top_k(self.feature_importances, self.k)
        return self
    def transform(self, X):
        return X[:, self.feature_indices_]
    

k = 6
top_k_feature_indices = indices_of_top_k(feature_importance, k)
np.array(attributes)[top_k_feature_indices]<jupyter_output><empty_output><jupyter_text>Randomized Search<jupyter_code>param_distribs = {
        'n_estimators': randint(low=1, high=200),
        'max_features': randint(low=1, high=8),
    }

rnd_search = RandomizedSearchCV(forest_reg, param_distributions=param_distribs,
                                n_iter=10, cv=5, scoring='neg_mean_squared_error', random_state=42)
rnd_search.fit(train_ds_x, train_ds_y)<jupyter_output><empty_output><jupyter_text>Evaulate Model<jupyter_code>final_model = grid_search.best_estimator_
test_x = test_ds.drop("median_house_value", axis = 1)
test_y = test_ds["median_house_value"].copy()

test_ds_x = full_pipeline.transform(test_ds)

predictions = final_model.predict(test_ds_x)
final_rmse = np.sqrt(metrics.mean_squared_error(test_y, predictions))
final_rmse
confidence = 0.95
squared_err = (test_y - predictions) ** 2
np.sqrt(t.interval(confidence, len(squared_err)-1, loc = squared_err.mean(),
                  scale = sem(squared_err)))
full_pipeline_with_predictor = pipeline.Pipeline([
    ('Preparation', full_pipeline),
    ('Random_forest_predictor', final_model)
])<jupyter_output><empty_output><jupyter_text>Save and Load Models<jupyter_code>my_model = full_pipeline_with_predictor

joblib.dump(my_model, 'my_RF_model.pkl')

my_model_loaded = joblib.load('my_RF_model.pkl')
n_prediction = my_model_loaded.predict(test_ds)
final_rmse_n = np.sqrt(metrics.mean_squared_error(test_y, n_prediction))
final_rmse_n<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.13727036190175826, 'fraction_numerical': 0.06788614681069148, 'mean_word_length': 5.721780604133546, 'pattern_counts': {'":': 0, '<': 60, '<?xml version=': 0, '>': 60, 'https://': 1, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '16', 'dir': 'jupyter-structured-clean-dedup', 'id': '527481', 'license': 'MIT', 'n_tokens_mistral': 5209, 'n_tokens_neox': 4572, 'n_words': 799, 'path': 'House_price_prediction/Housing_price_prediciton.ipynb', 'repo_name': 'Ziminli/ML_implementations', 'size': '12683'}