<jupyter_start><jupyter_text>Hands-on machine learning chapter 2: Housing projectDownload the data<jupyter_code>import os
import tarfile
import urllib

DOWNLOAD_ROOT = "https://raw.githubusercontent.com/ageron/handson-ml2/master/"
HOUSING_PATH = os.path.join("datasets", "housing")
HOUSING_URL = DOWNLOAD_ROOT + "datasets/housing/housing.tgz"

def fetch_housing_data(housing_url=HOUSING_URL, housing_path=HOUSING_PATH):
    os.makedirs(housing_path, exist_ok=True)
    tgz_path = os.path.join(housing_path, "housing.tgz")
    urllib.request.urlretrieve(housing_url, tgz_path)
    housing_tgz = tarfile.open(tgz_path)
    housing_tgz.extractall(path=housing_path)
    housing_tgz.close()
fetch_housing_data()<jupyter_output><empty_output><jupyter_text>Load the data using Pandas<jupyter_code>import pandas as pd

def load_housing_data(housing_path=HOUSING_PATH):
    csv_path = os.path.join(housing_path, "housing.csv")
    return pd.read_csv(csv_path)<jupyter_output><empty_output><jupyter_text>Top 5 rows of the data<jupyter_code>housing = load_housing_data()
housing.head()
housing.info()<jupyter_output><class 'pandas.core.frame.DataFrame'>
RangeIndex: 20640 entries, 0 to 20639
Data columns (total 10 columns):
 #   Column              Non-Null Count  Dtype  
---  ------              --------------  -----  
 0   longitude           20640 non-null  float64
 1   latitude            20640 non-null  float64
 2   housing_median_age  20640 non-null  float64
 3   total_rooms         20640 non-null  float64
 4   total_bedrooms      20433 non-null  float64
 5   population          20640 non-null  float64
 6   households          20640 non-null  float64
 7   median_income       20640 non-null  float64
 8   median_house_value  20640 non-null  float64
 9   ocean_proximity     20640 non-null  object 
dtypes: float64(9), object(1)
memory usage: 1.6+ MB<jupyter_text>ocean_proximity investigation: It's type object, so let's see what categories exist and how many districts belong    to each category:<jupyter_code>housing["ocean_proximity"].value_counts()
housing.describe()<jupyter_output><empty_output><jupyter_text>Show the data in a histogram<jupyter_code>#matplotlib inline  # only in a jupyter notebook commented out
import matplotlib.pyplot as plt
housing.hist(bins=50, figsize=(20,15))
plt.show()  # optional in jupyter notebook -> will display anyway when cell is executed<jupyter_output><empty_output><jupyter_text>create a testset Ratio: 20% of the data and set it aside<jupyter_code>import numpy as np

def split_train_test(data, test_ratio):
    # np.random.seed(42) # is there to create always the same shuffle of a testset
    shuffled_indices = np.random.permutation(len(data))
    test_set_size = int(len(data) * test_ratio)
    test_indices = shuffled_indices[:test_set_size]
    train_indices = shuffled_indices[test_set_size:]
    return data.iloc[train_indices], data.iloc[test_indices]
train_set, test_set = split_train_test(housing, 0.2)
len(train_set)
len(test_set)<jupyter_output><empty_output><jupyter_text>Using the method of creating hashes from identifiers to keep the testset consistent (otherwise the above solution would also work, however, it would break if we would update the dataset since we always have exact the same testset, no matter what updated data we get<jupyter_code>from zlib import crc32

def test_set_check(identifier, test_ratio):
    return crc32(np.int64(identifier)) & 0xffffffff < test_ratio * 2**32
def split_train_test_by_id(data, test_ratio, id_column):
    ids = data[id_column]
    in_test_set = ids.apply(lambda id_: test_set_check(id_, test_ratio))
    return data.loc[~in_test_set], data.loc[in_test_set]
housing_with_id = housing.reset_index()  # adds an 'index' column to the dataset
train_set, test_set = split_train_test_by_id(housing_with_id, 0.2, "index")<jupyter_output><empty_output><jupyter_text>If the row index is used as a unique identifier, one needs to make sure that if new data must be added to the end of the dataset. Also no row should ever get deleted. If this is not possible one needs to create an identifierout of the most stable features in the dataset.**in housing case** this is *longitude* and *latitude*<jupyter_code>housing_with_id["id"] = housing["longitude"] * 1000 + housing["latitude"]
train_set, test_set = split_train_test_by_id(housing_with_id, 0.2, "id")<jupyter_output><empty_output><jupyter_text>One can split the train and the testset also using **sklearn** python library<jupyter_code>from sklearn.model_selection import train_test_split

train_set, test_set = train_test_split(housing, test_size=0.2, random_state=42)<jupyter_output><empty_output><jupyter_text>The median income seems to be a very important attribute (page. 54) to predict median house prices.Create income category attribute with five categories<jupyter_code>housing["income_cat"] = pd.cut(housing["median_income"],
                              bins=[0, 1.5, 3.0, 4.5, 6.0, np.inf],
                              labels=[1, 2, 3, 4, 5])
housing["income_cat"].hist()<jupyter_output><empty_output><jupyter_text>Now one is ready to do stratified sampling based on the income category<jupyter_code>from sklearn.model_selection import StratifiedShuffleSplit

split = StratifiedShuffleSplit(n_splits=1, test_size=0.2, random_state=42)
for train_index, test_index in split.split(housing, housing["income_cat"]):
    strat_train_set = housing.loc[train_index]
    strat_test_set = housing.loc[test_index]
strat_test_set["income_cat"].value_counts() / len(strat_test_set)<jupyter_output><empty_output><jupyter_text>**Bring back data to normal:** Drop income_cat attribute<jupyter_code>for set_ in (strat_train_set, strat_test_set):
    set_.drop("income_cat", axis=1, inplace=True)<jupyter_output><empty_output><jupyter_text>Explore the data Create a copy<jupyter_code>housing = strat_train_set.copy()
housing.plot(kind="scatter", x="longitude", y="latitude", alpha=0.1)<jupyter_output><empty_output><jupyter_text>Let's look at the housing prices:<jupyter_code>housing.plot(kind="scatter", x="longitude", y="latitude", alpha=0.4,
            s=housing["population"]/100, label="population", figsize=(10,7),
             c="median_house_value", cmap=plt.get_cmap("jet"), colorbar=True)
plt.legend()<jupyter_output><empty_output><jupyter_text>looking for correlationsHow much does each feature correlate with the median house value<jupyter_code>corr_matrix = housing.corr()
corr_matrix["median_house_value"].sort_values(ascending=False)<jupyter_output><empty_output><jupyter_text>Correlation coefficient: ranges from -1 to 1. Where 0 means no correlation at all.From this example now we could take the following assumption: The median_house_value is highly correlating with the median_income in the way that "when the income is high, the median_house_value is also high"Another way to check for correlation between is the following `pandas scatter_matrix`<jupyter_code>from pandas.plotting import scatter_matrix

attributes = ["median_house_value", "median_income", "total_rooms", "housing_median_age"]
scatter_matrix(housing[attributes], figsize=(12,8))<jupyter_output><empty_output><jupyter_text>Most promising is the `median_income` towards `median_house_value` correlation:<jupyter_code>housing.plot(kind="scatter", x="median_income", y="median_house_value", alpha=0.1)<jupyter_output><empty_output><jupyter_text>Correlation is indeed very strong; you can clearly see the upwards trend when the income goes up, the median_house_value goes up as well  Experimenting with Attribute Combinations(page 61)Create new attributes `rooms_per_household`, `bedrooms_per_room` and `population_per_household`<jupyter_code>housing["rooms_per_household"] = housing["total_rooms"] / housing["households"]
housing["bedrooms_per_room"] = housing["total_bedrooms"] / housing["total_rooms"]
housing["population_per_housholds"] = housing["population"] / housing["households"]
corr_matrix = housing.corr()
corr_matrix["median_house_value"].sort_values(ascending=False)<jupyter_output><empty_output><jupyter_text>Prepare the Data for Machine Learning Algorithms(Page 62)<jupyter_code>housing = strat_train_set.drop("median_house_value", axis=1)
housing_labels = strat_train_set["median_house_value"].copy()<jupyter_output><empty_output><jupyter_text>the attribute `total_bedrooms` has some missing values, we need to fix this with one of the following options<jupyter_code>#housing.dropna(subset=["total_bedrooms"]) # Option1: get rid of the corresponding districts
housing.drop("total_bedrooms", axis=1) # Option2: get rid of the whole attribute
#median = housing["total_bedrooms"].median() # Option3: set the values to some value
#housing["total_bedrooms"].fillna(median, inplace=True)

from sklearn.impute import SimpleImputer

imputer = SimpleImputer(strategy="median")
# All numerical columns (so all columns except ocean_proximity)
housing_num = housing.drop("ocean_proximity", axis=1)
imputer.fit(housing_num)
imputer.statistics_
housing_num.median().values
# replace missing values with the learned medians of the trained imputer
X = imputer.transform(housing_num)
# back into pandas frame
housing_tr = pd.DataFrame(X, columns=housing_num.columns, index=housing_num.index)<jupyter_output><empty_output><jupyter_text>Handling Text and categorical attributes<jupyter_code>housing_cat = housing[["ocean_proximity"]]
housing_cat.head(10)
# the ocean_proximity is a categorical attribute. We turn those now into categories with scikit-learn
# Page 66
from sklearn.preprocessing import OrdinalEncoder
ordinal_encoder = OrdinalEncoder()
housing_cat_encoded = ordinal_encoder.fit_transform(housing_cat)
housing_cat_encoded[:10]
# list of categories 
ordinal_encoder.categories_
# Transform this into one-hot encoding strategy. Since enumerated categories let assume that two
# nearby values (0 and 1) are more similar than distant values. This is a wrong assumption here (<1H ocean is
# not closer to inland than near bay) -> solution: one-hot-encoding 
# Page 67
from sklearn.preprocessing import OneHotEncoder
cat_encoder = OneHotEncoder()
housing_cat_1hot = cat_encoder.fit_transform(housing_cat)
housing_cat_1hot
# If wanted to convert to a dense array (sparse array only stores positions of 1 instead of all zeros)
housing_cat_1hot.toarray()
cat_encoder.categories_
# Custom Transformers
# We create now these interesting fields that we looked at before: rooms_per_household, population_per_household
# and bedrooms_per_room (see step 29 (page 62))
from sklearn.base import BaseEstimator, TransformerMixin

rooms_ix, bedrooms_ix, population_ix, households_ix = 3, 4, 5, 6

class CombinedAttributesAdder(BaseEstimator, TransformerMixin):
    def __init__(self, add_bedrooms_per_room = True):  # no *args or **kwargs
        self.add_bedrooms_per_room = add_bedrooms_per_room
    def fit(self, X, y=None):
        return self  # nothing else to do
    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 = CombinedAttributesAdder(add_bedrooms_per_room=False)
housing_extra_attribs = attr_adder.transform(housing.values)<jupyter_output><empty_output><jupyter_text>Feature ScalingS. 70<jupyter_code># Transformation pipelines
from sklearn.pipeline import Pipeline
from sklearn.preprocessing import StandardScaler

num_pipeline = Pipeline([
    ('imputer', SimpleImputer(strategy="median")),
    ('attribs_adder', CombinedAttributesAdder()),
    ('std_scaler', StandardScaler())
])

housing_num_tr = num_pipeline.fit_transform(housing_num)
from sklearn.compose import ColumnTransformer

num_attribs = list(housing_num)
cat_attribs = ["ocean_proximity"]

full_pipeline = ColumnTransformer([
    ("num", num_pipeline, num_attribs),
    ("cat", OneHotEncoder(), cat_attribs),
])

housing_prepared = full_pipeline.fit_transform(housing)<jupyter_output><empty_output><jupyter_text>Select and train the model  Training and evaluating on the training setpage 72<jupyter_code># Linear regression
from sklearn.linear_model import LinearRegression

lin_reg = LinearRegression()
lin_reg.fit(housing_prepared, housing_labels)
# let's have a look at the data 
some_data = housing.iloc[:5]
some_labels = housing_labels.iloc[:5]
some_data_prepared = full_pipeline.transform(some_data)
print("Predictions: ", lin_reg.predict(some_data_prepared))
print("Labels:", list(some_labels))
from sklearn.metrics import mean_squared_error
housing_predictions = lin_reg.predict(housing_prepared)
lin_mse = mean_squared_error(housing_labels, housing_predictions)
lin_rmse = np.sqrt(lin_mse)
lin_rmse
# Decision Tree
from sklearn.tree import DecisionTreeRegressor

tree_reg = DecisionTreeRegressor()
tree_reg.fit(housing_prepared, housing_labels)
housing_predictions = tree_reg.predict(housing_prepared)
tree_mse = mean_squared_error(housing_labels, housing_predictions)
tree_rmse = np.sqrt(tree_mse)
tree_rmse<jupyter_output><empty_output><jupyter_text>Better evaluation using Cross-validationpage 73<jupyter_code># K-fold cross-validation feature (Scikit-learn)
# Splits randomly the training set into 10 distinct subsets (folds), then it trains and evaluates the Decision 
#  Tree model 10 times, picking a different fold for evaluation and every time and training on the other 
#  9 folds
# One real advantage is also the standard deviation which we only have if we train the model multiple times
#  (we would not have that when only training with one training set)
from sklearn.model_selection import cross_val_score
scores = cross_val_score(tree_reg, housing_prepared, housing_labels, 
                        scoring="neg_mean_squared_error", cv=10)
tree_rmse_scores = np.sqrt(-scores)
def display_scores(scores):
    print("Scores:", scores)
    print("Mean:", scores.mean())
    print("Standard deviation", scores.std())
display_scores(tree_rmse_scores)
# Cross validation for linear regression
lin_scores = cross_val_score(lin_reg, housing_prepared, housing_labels, 
                            scoring="neg_mean_squared_error", cv=10)
lin_rmse_scores = np.sqrt(-lin_scores)
display_scores(lin_rmse_scores)
# It seems that the decision tree model is overfitting very bad so that it's worse than the linear regression
# RandomForest Regressor
from sklearn.ensemble import RandomForestRegressor
forest_reg = RandomForestRegressor()
forest_reg.fit(housing_prepared, housing_labels)
housing_predictions = forest_reg.predict(housing_prepared)
forest_mse = mean_squared_error(housing_labels, housing_predictions)
forest_rmse = np.sqrt(forest_mse)
forest_rmse
forest_scores = cross_val_score(forest_reg, housing_prepared, housing_labels, 
                               scoring="neg_mean_squared_error", cv=10)
forest_rmse_scores = np.sqrt(-forest_scores)
display_scores(forest_rmse_scores)
# Dump the trained models:
import joblib

joblib.dump(lin_reg, "trained_models/linear_regression.pkl") # lin. regression
joblib.dump(tree_reg, "trained_models/decision_tree.pkl") # decision tree
joblib.dump(forest_reg, "trained_models/random_forest.pkl") # random forest<jupyter_output><empty_output><jupyter_text>Fine tune the modelp.76Fiddle with the hyperparameters using Scikit-learn's grid search. It tries out many combinations and cross-validates them<jupyter_code>from sklearn.model_selection import GridSearchCV

param_grid = [
    {'n_estimators': [3, 10, 30], 'max_features': [2, 4, 6, 8]},
    {'bootstrap': [False], 'n_estimators': [3, 10], 'max_features': [2, 3, 4]},
]

forest_reg = RandomForestRegressor()

grid_search = GridSearchCV(forest_reg, param_grid, cv=5,
                          scoring='neg_mean_squared_error',
                          return_train_score=True)

grid_search.fit(housing_prepared, housing_labels)
# directly get the best estimator of grid search
grid_search.best_estimator_
# evaluation scores
cvres = grid_search.cv_results_
for mean_score, params in zip(cvres["mean_test_score"], cvres["params"]):
    print(np.sqrt(-mean_score), params)
# analyze the best models and their errors
feature_importances = grid_search.best_estimator_.feature_importances_
feature_importances
# lets display the importance scores next to the corresponding attribute names
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_importances, attributes), reverse=True)<jupyter_output><empty_output><jupyter_text>Evaluate the system on the test setp. 79<jupyter_code>final_model = grid_search.best_estimator_

X_test = strat_test_set.drop("median_house_value", axis=1)
Y_test = strat_test_set["median_house_value"].copy()

X_test_prepared = full_pipeline.transform(X_test)

final_predictions = final_model.predict(X_test_prepared)

final_mse = mean_squared_error(Y_test, final_predictions)
final_rmse = np.sqrt(final_mse)  # => evaluates to 47730.2

# Confidence
from scipy import stats
confidence = 0.95
squared_errors = (final_predictions - Y_test) ** 2
np.sqrt(stats.t.interval(confidence, len(squared_errors) - 1,
                        loc=squared_errors.mean(),
                        scale=stats.sem(squared_errors)))<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.09224568574023614, 'fraction_numerical': 0.014702543142597639, 'mean_word_length': 6.109362389023406, 'pattern_counts': {'":': 0, '<': 107, '<?xml version=': 0, '>': 108, 'https://': 1, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 1, 'word_list_counts': {'cursed_substrings.json': 1, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '26', 'dir': 'jupyter-structured-clean-dedup', 'id': '115246', 'license': 'MIT', 'n_tokens_mistral': 5974, 'n_tokens_neox': 5602, 'n_words': 1368, 'path': 'Housing.ipynb', 'repo_name': 'pkhurt/hands-on-machine-learning', 'size': '17616'}