<jupyter_start><jupyter_text>Grabbing the Data<jupyter_code>def load_data(path="../datasets/housing"):
    csv_path = join(path, "housing.csv")
    return pd.read_csv(csv_path)
housing = load_data()
housing.head()<jupyter_output><empty_output><jupyter_text>Exploring the Data<jupyter_code>housing.info()
housing.hist(bins=50, figsize=(20, 15))
plt.show()
import hashlib

# Separate out the test set
# This will generate a new random test set every time though, so it's not optimal!
def split_data(data, ratio):
    # Mix up all of the data
    shuffled_indices = np.random.permutation(len(data))
    
    # Get the number of elements to put in the test set
    test_set_size = int(len(data) * ratio)
    
    # Separate out the test set and training set
    test_indices = shuffled_indices[:test_set_size]
    train_indices = shuffled_indices[test_set_size:]
    
    # Return the datasets
    return data.iloc[train_indices], data.iloc[test_indices]

# This helper function is used to check if data previously was in the training set
def test_set_check(identifier, ratio, hash):
    return hash(np.int64(identifier)).digest()[-1] < 256 * ratio

# This better split method will try to make sure that training data isn't used as testing data later
def split_by_id(data, ratio, id_col, hash=hashlib.md5):
    ids = data[id_col]
    in_set = ids.apply(lambda id_: test_set_check(id_, ratio, hash))
    return data.loc[~in_set], data.loc[in_set]
housing_with_id = housing.reset_index()
train_set, test_set = split_by_id(housing_with_id, 0.2, "index") # Use the better split function to get our data

housing_with_id['id'] = housing['longitude'] * 1000 + housing['latitude'] # Add 'index' column using
                                                                        # latitude and longitude
train_set, test_set = split_by_id(housing_with_id, 0.2, 'id')
print(len(train_set), "train +", len(test_set), "test")


#Alternatively, one can use the train_test_split function provided in sklearn.
#
#from sklearn.model_selection import train_test_split
#
#It basically works like our split_data function, but with extra params
#
#Random sampling is bad if we don't have a lot of data! (but you knew that from MATH 183)

# Let's add an income category attribute
housing['income_cat'] = np.ceil(housing['median_income'] / 1.5)
housing['income_cat'].where(housing['income_cat'] < 5, 5.0, inplace=True)

plt.hist(housing['income_cat'])
plt.show()
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)
# Let's remove income_cat so that we restore the data to it's original state
for set_ in (strat_train_set, strat_test_set):
    set_.drop("income_cat", axis=1, inplace=True)
housing = strat_train_set.copy()
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()

# It looks like california! (which makes sense, given what the data is...)
import matplotlib.image as mpimg
california_img=mpimg.imread('../images/end_to_end_project/california.png')
ax = housing.plot(kind="scatter", x="longitude", y="latitude", figsize=(10,7),
                       s=housing['population']/100, label="Population",
                       c="median_house_value", cmap=plt.get_cmap("jet"),
                       colorbar=False, alpha=0.4,
                      )
plt.imshow(california_img, extent=[-124.55, -113.80, 32.45, 42.05], alpha=0.5,
           cmap=plt.get_cmap("jet"))
plt.ylabel("Latitude", fontsize=14)
plt.xlabel("Longitude", fontsize=14)

prices = housing["median_house_value"]
tick_values = np.linspace(prices.min(), prices.max(), 11)
cbar = plt.colorbar()
cbar.ax.set_yticklabels(["$%dk"%(round(v/1000)) for v in tick_values], fontsize=14)
cbar.set_label('Median House Value', fontsize=16)

plt.legend(fontsize=16)
plt.show()
# Looking for Correlations
# Let's ge tthe standard correlation coefficient
corr_matrix = housing.corr()
corr_matrix['median_house_value'].sort_values(ascending=False)

# correlations are from [-1, 1]. Close the absolute value is to 1, the stronger the correlation.
# Works only for linear correlations
# Let's use the pandas scatter matrix to see the more strongly correlated values
# Using all 11 numerical attributes would give us 121 plots

attributes = ['median_house_value', 'median_income', 'total_rooms', 'housing_median_age']

pd.plotting.scatter_matrix(housing[attributes], figsize=(12, 8))
# Median house value seems strongly correlated to median income, so let's plot those more closely
housing.plot(kind='scatter', x='median_income', y='median_house_value', alpha=0.1)<jupyter_output><empty_output><jupyter_text>Prepping the Data<jupyter_code># First, revert back to a clean training set

# Drop the labels from our training set
housing = strat_train_set.drop('median_house_value', axis=1)

# Put the labels in a new Series
housing_labels = strat_train_set['median_house_value'].copy()<jupyter_output><empty_output><jupyter_text>Data CleaningUsually, there are three options to use to clean data1) Get rid of corresponding missing data2) Get rid of the entire attribute3) Set the values to some other value (zero, mean, median, etc)This can be done with Pandas methods easily1) `housing.dropna(subset=)`2) `housing.drop(, axis=1)`3) `housing[].fillna(, inplace=True)`*note: choosing option 3 means that the value to calculate should be calculated on the training set and usedon the training set, but should also be saved to replace the missing values in the test set as well*<jupyter_code># SKLearn has a nice class for handling missing values called Imputer. Steps to use it 
# 1) create an Imputer instance, specifying how you want to replace missing values
#   imputer needs numeric data only, so you'll need to create a copy of the data with only numeric values
# 2) Fit the instance to the data using .fit()
# 3) Use the "trained" imputer instance to transform the training set by replacing missing values by the learned
#   medians
# 4) The result is a numpy array with the transformed features, so create a dataframe out of it

from sklearn.preprocessing import Imputer

imputer = Imputer(strategy="median") #step 1
housing_nums = housing.drop('ocean_proximity', axis=1) #step 1
imputer.fit(housing_nums) #step 2
X = imputer.transform(housing_nums) #step 3
housing_tr = pd.DataFrame(X, columns=housing_nums.columns) #step 4
housing_tr.head()<jupyter_output><empty_output><jupyter_text>Handling Text and Categorical Attributes<jupyter_code>housing_cat = housing[['ocean_proximity']]
housing_cat.head(10)
# We need to turn this text into numbers for our algorithms to actually use
# For this, we can use the .factorize() method from Pandas which will map each
# category to a different integer
from future_encoders import OrdinalEncoder

ord_enc = OrdinalEncoder()
housing_cat_encoded = ord_enc.fit_transform(housing_cat)
housing_cat_encoded[:10]

# A potential issue with this is that ML algos will assume that two nearby values are more similar
# than two distant values. To fix this, a common solution is to create one binary attribute per
# category. this is called one-hot encoding because only one attribute will be equal to 1 (hot) while
# the other will be 0 (cold)
# To do the aforementioned encoding, we can use SkLearn's OneHotEncoder
from future_encoders import OneHotEncoder
cat_encoder = OneHotEncoder()
housing_cat_1hot = cat_encoder.fit_transform(housing_cat)
housing_cat_1hot<jupyter_output><empty_output><jupyter_text>Writing Custom TransformersThis may be necessary to do more complex transformations, so creating a class that works semlessly with SkLearn's functionalities is imperative. All you need to do this is to create a class and implement the `fit()` (returning `self`), `transform()`, and `fit_transform()` methods. The last one can be gotten for free by extending from the `TransformerMixin` class. If you also extend from `BaseEstimator` (and don't use `*args` and `**kargs` in your constructor) you can get two extra methods (`get_params()` and `set_params()`) that will be useful for hyperparameter tuning. The following is a sample class as an example.```from sklearn.base import BaseEstimator, TransformerMixinrooms_ix, bedrooms_ix, population_ix, household_ix = 3, 4, 5, 6class CombinedAttributesAdder(BaseEstimator, 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, y=None):        rooms_per_household = X[:, rooms_ix] / X[:, household_ix]        population_per_household = X[:, population_ix] / X[:, household_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] ```  This class can be used like so:  ``` attr_adder = CombinedAttributesAdded(False) housing_extra_attribs = attr_adder.transform(housing.values) ```<jupyter_code>from sklearn.base import BaseEstimator, TransformerMixin

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

class CombinedAttributesAdder(BaseEstimator, 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, y=None):
        rooms_per_household = X[:, rooms_ix] / X[:, household_ix]
        population_per_household = X[:, population_ix] / X[:, household_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)
housing_extra_attribs = pd.DataFrame(
    housing_extra_attribs,
    columns=list(housing.columns)+["rooms_per_household", "population_per_household"])
housing_extra_attribs.head()<jupyter_output><empty_output><jupyter_text>Feature ScalingIt's important to make sure all of the data has the same scale for the algos to work properly. There are two ways to do this: _**min-max scaling (normalization)**_ and _**standardization**_The former is simple: just shift all values and scale them so that the range from [0, 1]. This is done by subtracting the min value and dividing by max - min. SKLearn gives a Transformer called `MinMaxScaler` for this (with a `feature_range` hyperparameter that lets you change the range if you don't want 0-1)The latter subtracts the mean value (so all standardized values have a 0 mean) then divides by the variance so that the resulting distribution has unit variance. It doesn't bound to a specific range, which may be a problem for some algos (like neural nets). However, it's also much less affected by outliers in the data. SKLearn has a Transformer called the `StandardScaler` for this.*note: make sure to fit scalers **only** to the training data! After this, you can use them to transform both the training set and testing set, as well as new data* Transformation PipelinesMany of the steps for Data Transformation need to occur in the right order. SKLearn has a `Pipeline` class which helps with such sequences. For the numerical attributes of our data, we can do this<jupyter_code># Build a pipeline for the numerical attributes of our data
from sklearn.pipeline import Pipeline
from sklearn.preprocessing import StandardScaler

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

housing_num_tr = num_pipeline.fit_transform(housing_nums)
housing_num_tr
# Let's now transform our DataFrame into a Numpy array for use

from sklearn.base import BaseEstimator, TransformerMixin

class DataFrameSelector(BaseEstimator, TransformerMixin):
    def __init__(self, attribute_names):
        self.attribute_names = attribute_names
    def fit(self, X, y=None):
        return self
    def transform(self, X):
        return X[self.attribute_names].values
# Using this class, let's make a pipleine that will take the dataframe and handle the numeric values
num_attribs = list(housing_nums)
cat_attribs = ['ocean_proximity']

num_pipeline = Pipeline([
    ('selector', DataFrameSelector(num_attribs)),
    ('imputer', Imputer(strategy='median')),
    ('attribs_adder', CombinedAttributesAdder()),
    ('std_scaler', StandardScaler())
])

cat_pipeline = Pipeline([
    ('selector', DataFrameSelector(cat_attribs)),
    ('cat_encoder', OneHotEncoder(sparse=False))
])
# Link our pipelines up
from sklearn.pipeline import FeatureUnion

full_pipeline = FeatureUnion(transformer_list=[
    ('num_pipeline', num_pipeline),
    ('cat_pipeline', cat_pipeline)
])
housing_prepared = full_pipeline.fit_transform(housing)
housing_prepared
housing_prepared.shape<jupyter_output><empty_output><jupyter_text>Selecting and Training a Model Training and Evaluating on the Training Set<jupyter_code># Let's first do a LinearRegression model

from sklearn.linear_model import LinearRegression

lin_reg = LinearRegression()
lin_reg.fit(housing_prepared, housing_labels)
# Let's try it out
some_data = housing.iloc[:5]
some_labels = housing_labels.iloc[:5]
some_data_prepped = full_pipeline.transform(some_data)
print("Predictions: {}".format(lin_reg.predict(some_data_prepped)))
print("Labels: {}".format(list(some_labels)))
# Issues
# The models is really, really off. That's not good at all. Let's check the RMSE
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)
print(lin_rmse)

# Typical prediction error of $69K is not very good at all...
# This model underfits the training data. This is remedied by
# 1) selecting a more powerful model
# 2) using better features
# 3) reduce constraints on the model (not available for this model)
# Let's now try a DecisionTreeRegression model
from sklearn.tree import DecisionTreeRegressor

tree_reg = DecisionTreeRegressor()
tree_reg.fit(housing_prepared, housing_labels)
# Let's try this one out
housing_predictions = tree_reg.predict(housing_prepared)
tree_mse = mean_squared_error(housing_labels, housing_predictions)
tree_rmse = np.sqrt(tree_mse)
tree_rmse

# This model overfits the data. There's no way that we have 0 error at all<jupyter_output><empty_output><jupyter_text>Better Evaluation Using Cross-Validation<jupyter_code># Let's cross validate to better evaluate our decision tree model. We could do this manually by splitting
# the training set and validation set, then training against the smaller sets, but we can use SKLearn's built in
# classes instead

# The following code will do K-fold cross validation
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)

# The reason for the negative scores is that sklearn expects a utility function (greater == better) instead
# of a cost function (lower == better), so scoring is the opposite of MSE. This is why we negate scores
def display_scores(scores):
    print("Scores: {}".format(scores))
    print("Mean: {}".format(scores.mean()))
    print("Standard Deviation: {}".format(scores.std()))
    
display_scores(tree_rmse_scores)

# Our average error is now $71K instead of $69K, so this model is actually worse!
# Let's do the same analysis with the LinearRegression model, just to be sure
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)
# Let's try another model: RandomForestRegression
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)

# Better, but our model still overfits since the RMSE on the train data of $22K is substantially larger than
# the value in the test data ($52K) though the $52K is the best one so far

# From here, our goal would be to try out many models from other categories of ML algos without spending too much
# time tweaking hyperparams. The goal is to get 2-5 promising models for use
# Try SupportVectorRegression, now
from sklearn.svm import SVR

svm_reg = SVR(kernel="linear")
svm_reg.fit(housing_prepared, housing_labels)
housing_predictions = svm_reg.predict(housing_prepared)
svm_mse = mean_squared_error(housing_labels, housing_predictions)
svm_rmse = np.sqrt(svm_mse)
svm_rmse

# Overall, pretty bad. Severly underfits the data<jupyter_output><empty_output><jupyter_text>*note: to save models we use, use the `pickle` module or `sklearn.externals.joblib` to make it easy.*```from sklearn.externals import joblibjoblib.dump(my_model, "my_model.pkl")```*then later on...*```my_model_loaded = joblib.load("my_model.pkl")```  Fine Tuning Grid Search<jupyter_code># This method of tuning is to use SKLearn's GridSearchCV to try out combinations so that you don't
# manually search through/try a bunch of them

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')
grid_search.fit(housing_prepared, housing_labels)

# This code will explore 18 combinations of RandomForest hyperparams and will train each model 5 times
# due to the 5-fold cross validation. 90 rounds of training will be performed
# GridSearch is great for relatively few combinations, but when the hyperparam search space is large,
# it's usually better to use RandomizedSearchCV instead
print(grid_search.best_params_, "\n", 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)<jupyter_output>64468.94797504932 {'max_features': 2, 'n_estimators': 3}
55952.070347143825 {'max_features': 2, 'n_estimators': 10}
52897.610697221826 {'max_features': 2, 'n_estimators': 30}
60439.20424121991 {'max_features': 4, 'n_estimators': 3}
52592.23367500706 {'max_features': 4, 'n_estimators': 10}
50504.98371680477 {'max_features': 4, 'n_estimators': 30}
58985.9868955743 {'max_features': 6, 'n_estimators': 3}
52299.01092591614 {'max_features': 6, 'n_estimators': 10}
49866.943443655495 {'max_features': 6, 'n_estimators': 30}
58257.10306852905 {'max_features': 8, 'n_estimators': 3}
52216.5063402653 {'max_features': 8, 'n_estimators': 10}
50237.442879438575 {'max_features': 8, 'n_estimators': 30}
62857.8060811324 {'bootstrap': False, 'max_features': 2, 'n_estimators': 3}
53941.874477775455 {'bootstrap': False, 'max_features': 2, 'n_estimators': 10}
59307.66510381341 {'bootstrap': False, 'max_features': 3, 'n_estimators': 3}
52518.36445081219 {'bootstrap': False, 'max_features': 3, 'n_estimators': [...]<jupyter_text>RandomizedSearchCV<jupyter_code># This option has the benefits of being able to explore a lot more hyperparams
# You have more control of computing power used by setting the number of iterations<jupyter_output><empty_output><jupyter_text>Ensemble Methods<jupyter_code># Combining models that perform the best will often perform better than the best individual model
# Think of how RandomForests was better than regular Decision Trees
# Let's check the relative importances of each attribute
feature_importances = grid_search.best_estimator_.feature_importances_
feature_importances
extra_attribs = ['rooms_per_hhold', 'pop_per_hhold', 'bedrooms_per_room']
cat_encoder = cat_pipeline.named_steps['cat_encoder']
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>Knowing this info, we could drop the less useful categories (like how only one `ocean_proximity` category is useful so the rest may be worth dropping  Evaluate System on the Test Set<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)
final_rmse
# Compute a 95% confidence interval with t-scores instead of z-scores (t distribution vs normal distribution)
from scipy import stats

confidence = 0.95
squared_errors = (final_predictions - y_test) ** 2
mean = squared_errors.mean()
m = len(squared_errors)

np.sqrt(stats.t.interval(confidence, m - 1,
                        loc=np.mean(squared_errors),
                        scale=stats.sem(squared_errors)))<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.09380219383379226, 'fraction_numerical': 0.02234460343971471, 'mean_word_length': 5.523557126030624, 'pattern_counts': {'":': 0, '<': 54, '<?xml version=': 0, '>': 52, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 3, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '13', 'dir': 'jupyter-structured-clean-dedup', 'id': '87220', 'license': 'Apache-2.0', 'n_tokens_mistral': 7423, 'n_tokens_neox': 6833, 'n_words': 2162, 'path': 'exercises/Chapter 2 - End-to-End Machine Learning Project.ipynb', 'repo_name': 'nate-browne/Intro-to-ML', 'size': '22153'}