此文章参考《hands on machine learning with sklearn and tensorflow》一书， 以加州房价数据集， 并建立加州房价模型的过程， 以此来记录一个完整机器学习项目的完整流程。
与原书有点出入的地方：原书是先进行数据集划分，再进行数据清洗操作。而该文章中则是先进行数据清洗，再进行数据集划分。其理由如下：
1、无论训练集和测试集都需要进行数据清洗操作才能更好的被机器学习算法处理；
2、在大部分情况下，这两组操作的顺序不同，没有什么问题。但是如果像在自然语言处理中，需要删除频率低的词时，先划分数据集的话，会造成训练集和测试集的数据分布不一致。
综上，所以一般习惯性是先对整个数据集进行数据清洗，然后再进行随机数据集划分。

1、项目概述

利用加州普查数据，建立一个加州房价模型。这个数据包含每个街区组的人口、收入中位数、房价中位数等指标。你的模型要利用这个数据进行学习，然后根据新的数据的指标，预测其对应街区的房价中位数。

1.1 划定问题

通常一个完整的机器学习工程由多个机器学习模型组成，其中一个模型的输入可能是另一个模型的输出，故第一步我们需要先弄清楚该模型的输入和输出分别是什么。其次，需要划定问题：监控或非监督，还是强化学习？这是个分类问题、回归问题还是其它？要使用批量学习还是线上学习？

在该项目中，模型的输入是一个加州普查得到的关于房价的数据信息，输出是各个街区组的房价中位数，则该模型的学习方式属于监督学习，且为多变量回归算法。

2、获取数据

数据下载代码如下：

import os
import tarfile
from six.moves import urllib

DOWNLOAD_ROOT = "https://raw.githubusercontent.com/ageron/handson-ml/master/"
HOUSING_PATH = "datasets/housing"
HOUSING_URL = DOWNLOAD_ROOT + HOUSING_PATH + "/housing.tgz"

def fetch_housing_data(housing_url=HOUSING_URL, ousing_path=HOUSING_PATH):

    if not os.path.isdir(housing_path):
      os.makedirs(housing_path)

    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()

数据集加载代码如下：

import pandas as pd
import os

def load_housing_data(housing_path=HOUSING_PATH):
    csv_path = os.path.join(housing_path, "housing.csv")
    return pd.read_csv(csv_path)

3、数据探索和可视化，发现规律

这部门内容主要基于 Pandas 库来实现。

3.1 快速查看数据结构

housing = load_housing_data()
# DataFrame 对象自带的 head() 函数可以查看数据集的前 5 行数据
housing.head()

前5行数据.png

每一行表示一个街区，一个样本数据。共有10个属性：经度、维度、房屋年龄中位数、总房间数、总卧室数、人口数、家庭数、收入中位数、房屋价值中位数、离大海距离。

# info() 方法可以快速查看数据的描述，特别是总行数、每个属性的类型和非空值的数量
housing.info()

info.png

从上图可以看到数据集共有20640个实例，但总卧室数只有20433个非空值，意味着有207个街区缺少这个值。

# 查询 ocean_proximity 属性的取值有多少个，及每个取值各有多少个样本数据
housing["ocean_proximity"].value_counts()

ocean_proximity.png

# describe() 方法展示了数值属性的概括
housing.describe()

数据属性概述.png

上图所示，describe() 方法统计了各个属性值的数量、均值、最小值和最大值等内容

另一种快速了解数据类型的方法是画出每个数值属性的柱状图。

%matplotlib inline  # only in a Jupyter notebook
import matplotlib.pyplot as plt
housing.hist(bins=50, figsize=(20, 15))
plt.show()

各个属性的柱状图.png

3.2 地理数据可视化

这个数据集比较特殊，里面包含了经纬度，我们可以创建一个所有街区的散点图来数据可视化：

# 依赖于 matplotlib
housing.plot(kind="scatter", x="longitude", y="latitude")

bad_visualization_plot.png

# 第二个方式展示可以更好展示数据密集性
housing.plot(kind="scatter", x="longitude", y="latitude", alpha=0.1)

better_visualization_plot.png

接下来展示房价、位置和人口密度的联系关系：

housing.plot(kind="scatter", x="longitude", y="latitude", alpha=0.4, 
                    s=housing["population"]/100, label="population", 
                    c="median_house_value", cmap=plt.get_cmap("jet"), 
                    colorbar=True, sharex=False)
plt.legend()

housing_prices_scatterplot.png

3.3 查找关联

# 因为数据集并不是非常大，可以很容易地使用 corr() 方法计算出每对属性间的标准相关系数（也称作皮尔逊相关系数）
corr_matrix = housing.corr()
corr_matrix["median_house_value"].sort_values(ascending=False)

median_house_value 的相关系数.png

另一种检测属性间相关系数的方法是使用 Pandas 的 scatter_matrix 函数，能画出每个数值属性对每个其它数值属性的图。此处指关注几个和房价中位数最有可能相关的属性：

from pandas.tools.plotting import scatter_matrix
attributes = ["median_house_value", "median_income", "total_rooms", "housing_median_age"]
scatter_matrix( housing[attributes], figsize=(12,8) )

scatter_matrix 示意图.png

3.4 属性组合试验

根据已有的属性数据计算或组合形成新的属性数据：
结论：与总房间数和卧室数相比，新的 bedrooms_per_room 属性与房价中位数的关联更强。

housing["rooms_per_household"] = housing["total_rooms"]/housing["households"]
housing["bedrooms_per_room"] = housing["total_bedrooms"]/housing["total_rooms"]
housing["population_per_household"]=housing["population"]/housing["households"]

corr_matrix = housing.corr()
corr_matrix["median_house_value"].sort_values(ascending=False)

median_house_value 的新的相关系数.png

4、为机器学习算法准备数据

编写一些函数来实现数据处理，不用手工来做的理由如下：

• 函数可以让你在任何数据集上方便的进行重复数据转换；
• 你能慢慢建立一个转换函数库，可以在未来的项目中复用；
• 在将数据传给算法之前，你可以在实时系统中使用这些函数；
• 这可以让你方便地尝试多种数据转换，查看哪些转换方法结合起来效果最好。

4.1 数据清洗

4.1.1 处理数值属性

大多机器学习算法不能处理缺失的特征，因此先创建一些函数来处理特征缺失的问题，此处采用普通的办法：

• 去掉对应的街区：

housing.dropna(subset=["total_bedrooms"])

• 去掉整个属性：

housing.drop("total_bedrooms", axis=1)

• 进行赋值（0， 平均值，中位数等）：

median = housing["total_bedrooms"].median()
housing["total_bedrooms"].fillna(median)

采用 Scikit-Learn 提供了一个方便的类来处理缺失值：Imputer。具体示例如下：

• 指定用某属性的中位数来替换该属性所有的缺失值

# 这个方法只适用数值属性，不适用于文本属性
from sklearn.preprocessing import Imputer
imputer = Imputer(strategy="median")
# 删除文本属性列，只对数值属性列进行转换
housing_num = housing.drop("ocean_proximity", axis=1)
# 计算中位数
imputer.fit(housing_num)
# imputer 计算出了每个属性的中位数，并将结果保存在了实例变量 statistics_ 中
imputer.statistics_
# 使用这个“训练过的” imputer 来对训练集进行转换，将缺失值替换为中位数：
X = imputer.transform(housing_num)
# 上面装换后特征的普通的 Numpy 数组。如果你想将其放回到 Pandas DataFrame 中
housing_tr = pd.DataFrame(X, columns=housing_num.columns)

4.2 处理文本和类别属性

sklearn 提供了两个转换器：LabelEncoder 和 CategoricalEncoder:

• LabelEncoder : 该类应用于标签列的转换，也适用于单列的文本属性列；
• CategoricalEncoder : 该类应用于多列的文本属性列，但目前sklearn库里仍不包含

备注：此处附上CategoricalEncoder 类的源码：

from sklearn.base import BaseEstimator, TransformerMixin
from sklearn.utils import check_array
from sklearn.preprocessing import LabelEncoder
from scipy import sparse

class CategoricalEncoder(BaseEstimator, TransformerMixin):
    
    def __init__(self, encoding='onehot', categories='auto', dtype=np.float64, handle_unknown='error'):
        self.encoding = encoding
        self.categories = categories
        self.dtype = dtype
        self.handle_unknown = handle_unknown

    def fit(self, X, y=None):
        if self.encoding not in ['onehot', 'onehot-dense', 'ordinal']:
            template = ("encoding should be either 'onehot', 'onehot-dense'  or 'ordinal', got %s")
            raise ValueError(template % self.handle_unknown)
        if self.handle_unknown not in ['error', 'ignore']:
            template = ("handle_unknown should be either 'error' or 'ignore', got %s")
            raise ValueError(template % self.handle_unknown)
        if self.encoding == 'ordinal' and self.handle_unknown == 'ignore':
            raise ValueError("handle_unknown='ignore' is not supported for encoding='ordinal'")
        X = check_array(X, dtype=np.object, accept_sparse='csc', copy=True)
        n_samples, n_features = X.shape

        self._label_encoders_ = [LabelEncoder() for _ in range(n_features)]
        for i in range(n_features):
            le = self._label_encoders_[i]
            Xi = X[:, i]
            if self.categories == 'auto':
                le.fit(Xi)
            else:
                valid_mask = np.in1d(Xi, self.categories[i])
                if not np.all(valid_mask):
                    if self.handle_unknown == 'error':
                        diff = np.unique(Xi[~valid_mask])
                        msg = ("Found unknown categories {0} in column {1} during fit".format(diff, i))
                        raise ValueError(msg)
                le.classes_ = np.array(np.sort(self.categories[i]))
        self.categories_ = [le.classes_ for le in self._label_encoders_]

        return self

    def transform(self, X):
        X = check_array(X, accept_sparse='csc', dtype=np.object, copy=True)
        n_samples, n_features = X.shape
        X_int = np.zeros_like(X, dtype=np.int)
        X_mask = np.ones_like(X, dtype=np.bool)

        for i in range(n_features):
            valid_mask = np.in1d(X[:, i], self.categories_[i])
            
            if not np.all(valid_mask):
                if self.handle_unknown == 'error':
                    diff = np.unique(X[~valid_mask, i])
                    msg = ("Found unknown categories {0} in column {1} during transform".format(diff, i))
                    raise ValueError(msg)
                else:
                    X_mask[:, i] = valid_mask
                    X[:, i][~valid_mask] = self.categories_[i][0]

            X_int[:, i] = self._label_encoders_[i].transform(X[:, i])

        if self.encoding == 'ordinal':
            return X_int.astype(self.dtype, copy=False)

        mask = X_mask.ravel()
        n_values = [cats.shape[0] for cats in self.categories_]
        n_values = np.array([0] + n_values)
        indices = np.cumsum(n_values)

        column_indices = (X_int + indices[:-1]).ravel()[mask]
        row_indices = np.repeat(np.arange(n_samples, dtype=np.int32), n_features)[mask]
        data = np.ones(n_samples * n_features)[mask]

        out = sparse.csc_matrix((data, (row_indices, column_indices)), shape=(n_samples, indices[-1]), dtype=self.dtype).tocsr()

        if self.encoding == 'onehot-dense':
            return out.toarray()
        else:
            return out

LabelEncoder 使用示例：

from sklearn.preprocessing import LabelEncoder

encoder = LabelEncoder()
housing_cat = housing["ocean_proximity"]
housing_cat_encoded = encoder.fit_transform(housing_cat)

sklearn 提供了一个编码器 OneHotEncoder ，用于将整数分类值转变为独热向量。
备注：fir_transform() 用于 2D 数组，而 housing_cat_encoded 是一个 1D 数组

from sklearn.preprocessing import OneHotEncoder

encoder = OneHotEncoder()
# 输出结果是一个 SciPy 稀疏矩阵，而不是 NumPy 数组。
housing_cat_1hot = encoder.fit_transform(housing_cat_encoded.reshape(-1, 1))

onehot.png

4.3 自定义转换器

示例：属性组合的自定义转换器

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): # no *args or **kargs
        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, y=None):
        rooms_per_houshold = 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)

4.4 特征缩放

当输入的各个数值属性量度不同时，机器学习算法的性能都不会好，进行特征缩放有两个方法，分别是归一化和标准化：

• 归一化：减去最小值，然后再除以最大值与最小值的差值。sklearn 提供了一个转换器 MinMaxScaler 来实现这个功能
• 标准化：减去平均值，然后除以方差，使得到的分布具有单位方差。sklearn 提供了一个转换器 StandardScaler 来实现这个功能

4.5 转换流水线

sklearn 提供了类 Pipline，来进行这一系列的转换。
    Pipline 构造器需要一个定义步骤顺序的名字/估计器对的列表。除了最后一个估计器，其余都要是转换器（即，它们都要有 fit_transform() 方法）
    当调用流水线的 fit() 方法，就会对所有转换器顺序调用 fit_transform() 方法，将每次调用的输出作为参数传递给下一个调用，一直到最后一个估计器，它只执行 fit() 方法。
如下是一个数值属性的小流水线：

from sklearn.pipline import Pipline
from sklearn.preprocessing import StandardScaler
from sklearn.base import BaseEstimator, TransformerMixin

num_pipline = Pipline([
                      ('imputer', Imputer(strategy="median")),
                      ('attribs_adder', CombinedAttributesAdder()),
                      ('std_scaler', StandardScaler()),
                      ])
housing_num_tr = num_pipline.fit_transform(housing_num)

转换流水线结果截图.png

sklearn 提供了一个类 FeatureUnion 实现将各个流水线的输出结果合并成一个结果。如将数值属性的转换结果和文本属性的转换结果合并为一个结果，示例代码如下：

from sklearn.pipeline import FeatureUnion

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

# 用于选择数组属性 or 文本属性 的特征列
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

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

cat_pipeline = Pipeline([
                        ('selector', DataFrameSelector(cat_attribs)),
                        ('label_binarizer', LabelBinarizer()),
                      ])

full_pipeline = FeatureUnion(transformer_list = [
                        ("num_pipeline", num_pipeline),
                        ("cat_pipeline", cat_pipeline),
                      ])
housing_prepared = full_pipeline.fit_transform(housing)
housing_prepared

housing_prepared.png

4.6 划分数据集

传统机器学习项目一般将数据集划分为训练集和测试集，其中训练集用模型参数的学习，测试集则是用于验证模型的优劣。另外在深度学习中，还会额外划分验证集，验证集主要使用于调整模型的超参数。【超参数：属于人工调整参数，不是训练学习所得】

    关于数据集划分比例建议如下：如果数据集样本数量在万或以下级别，则一般按80%：20%进行划分；如果数据集样本数量在百万级别及以上，则随机抽取千级别数量的样本组成测试集即可。

sklearn 提供了 train_test_split 方法进行数据集划分方法，如下为示例：

from sklearn.model_selection import train_test_split
# 参数说明：
# housing : 需要划分的数据集
# test_size, train_size : 训练集和测试集划分的比例
# random_state : 随机生成器种子，设定后，重新运行代码划分结果一样的
train_set, test_set = train_test_split(housing, test_size=0.2, random_state=42)

分层抽样：如性别比例、收入分类等，其中 sklearn 提供了 StratifiedShuffleSplit 类：
如下示例：现将收入分类，再进行分层抽样：

housing["income_cat"] = np.ceil(housing["median_income"] / 1.5)
housing["income_cat"].where(housing["income_cat"] < 5, 5.0, inplace=True)

from sklearn.model_selection import StratifiedShuffleSplit

split = StratifiedShuffleSplit(n_split=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]

housing["income_cat"].value_counts() / len(housing)

income_cat.png

5、选择模型，进行训练

完成了数据清洗，并划分好数据集，接下来我们就可以选择合适的模型，并采用训练集的数据对其进行训练学习，并测试集进行验证测试。最终挑选合适的模型。

前面项目概述中划定问题中，这是一个监督学习中的回归问题，故我们可以选择的模型有：线性模型、CART和随机森林等。因为此处只是梳理机器学习项目的流程，所以我们这里选择线性模型进行训练：

from sklearn.linear_model import LinearRegression

lin_reg = LinearRegression()
lin_reg.fit(housing_prepared, housing_labels)
# 接下来用一些数据来验证训练的结果
some_data = housing.iloc[:5]
some_labels = housing_labels.iloc[:5]
some_data_prepared = full_pipeline.transform(some_data)

print("Predictions:\t", lin_reg.predict(some_data_prepared))
# Predictions: [ 303104. 44800. 308928. 294208. 368704.]
print("Labels:\t\t", list(some_labels))
# Labels: [359400.0, 69700.0, 302100.0, 301300.0, 351900.0]

上述代码是采用了sklearn提供的线性模型对数据集进行训练学习，如果上述模型不满足我们的需求，我们也可以自行实现一个线性模型，构建对应于的优化目标，即代价函数，并采用梯度下降优化方法来进行参数的优化学习，如果数据量小的情况下，我们也可以直接使用最小二乘法来直接求解参数的最优解。

下面我们使用 sklearn 的 mean_squared_error 函数，用全部训练集来计算下这个回归模型的 RMSE：

from sklearn.metrics import mean_squared_error

housing_predictions = lin_reg.predict(housing_prepared)
lin_mes = mean_squared_error(housing_labels, housing_predictions)
lin_rmse = np.sqrt(lin_mse)
lin_rmse
# 68628.413493824875

从上面的结果可以看出这是一个模型欠拟合训练数据的模型。基本可以确定线性模型无法很好拟合该数据集的分布。此时，最好的做法就是更换模型。

另一种更好衡量模型拟合程度的方法：交叉验证，如常用的 K 折交叉验证。
它随机地将训练集分成十个不同的子集，成为“折”，然后训练评估决策模型 K 次，每次选一个不用的折来做评估，用其它 K-1 个来做训练。结果是一个包含 10 个评分的数组。示例代码如下：

from sklearn.model_selection import cross_val_score

scores = cross_val_score(lin_reg, housing_prepared, housing_labels, 
                          scoring="neg_mean_squared_error", cv=10)

rmse_scores = np.sqrt(-scores)

6、微调模型

模型超参数微调的方法：

• 手工调整：耗时长，且非常依赖人的经验
• 网格搜索：由人设定每个超参数中的取值集合，然后由程序自行组合来运行模型，最后挑选最好的参数组合，代码示例如下：

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)
# 最佳的参数组合
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)

• 随机搜索

当超参数的搜索空间相对较少时，网格搜索还可以，但是如果当超参数的搜索空间很大时，最好使用随机搜索。即 sklearn 中的方法：RandomizedSearchCV。随机搜索有两个优点：
1、如果你让随机搜索运行，比如1000次，它会探索每个超参数的1000个不同的值；
2、你可以方便地通过设定搜索次数，控制超参数搜索的计算量。

• 集成方法

将表现最好的模型组合起来。组合之后的性能通常要比单独的模型要好。其中随机森林就是决策树的一个集成表示方式。

7、给出解决方案

1、分析最佳模型和它们的误差：通过分析最佳模型及对应的误差，可以更好进行参数调整及其它调整。如随机森林可以判断各个属性值的重要性，这样可以在数据清洗阶段去掉什么信息的特征等。
2、用测试集评估系统：采用测试集，对最佳模型进行评测，通过模型在测试集上的表现，来判定模型是否过拟合，以及误差和性能是否满足目标等。

8、部署、监控、维护系统

这一步就是将训练好且误差在可接受的范围内的模型部署上线，供其它生产系统使用等。其中就是设计的监控、测试等代码的配套设计和编写。

参考资料：《hands on machine learning with sklearn and tensorflow》