Python-数据科学秘籍-五-

Python 数据科学秘籍（五）

原文：annas-archive.org/md5/a4f348a4e11e27ea41410c793e63daff

译者：飞龙

协议：CC BY-NC-SA 4.0

第九章：生长树

本章我们将涵盖以下食谱：

• 从树到森林——随机森林

• 生长极度随机化的树

• 生长旋转森林

介绍

在本章中，我们将看到更多基于树的算法的袋装方法。由于它们对噪声的鲁棒性以及对各种问题的普适性，它们在数据科学社区中非常受欢迎。

大多数这些方法的名声在于它们相比其他方法能够在没有任何数据准备的情况下获得非常好的结果，而且它们可以作为黑盒工具交给软件工程师使用。

除了前文提到的过高的要求外，还有一些其他优点。

从设计上看，袋装法非常适合并行化。因此，这些方法可以轻松应用于集群环境中的大规模数据集。

决策树算法在树的每一层将输入数据划分为不同的区域。因此，它们执行了隐式的特征选择。特征选择是构建良好模型中的一个重要任务。通过提供隐式特征选择，决策树相较于其他技术处于有利位置。因此，带有决策树的袋装法具备这一优势。

决策树几乎不需要数据准备。例如，考虑属性的缩放。属性的缩放对决策树的结构没有影响。此外，缺失值不会影响决策树。异常值对决策树的影响也很小。

在我们之前的一些食谱中，我们使用了多项式特征，仅保留了交互项。通过集成树方法，这些交互关系得到了处理。我们无需进行显式的特征转换来适应特征交互。

基于线性回归的模型在输入数据中存在非线性关系时会失败。当我们解释核主成分分析（Kernel PCA）食谱时，我们看到过这种效果。基于树的算法不受数据中非线性关系的影响。

对于树基方法的主要投诉之一是树的剪枝困难，容易导致过拟合。大树往往也会拟合底层数据中的噪声，从而导致低偏差和高方差。然而，当我们生长大量树木，并且最终预测是所有树的输出的平均值时，就能避免方差问题。

本章我们将介绍三种基于树的集成方法。

我们的第一个食谱是实现随机森林用于分类问题。Leo Breiman 是这一算法的发明者。随机森林是一种集成技术，通过内部使用大量的树来构建模型，用于解决回归或分类问题。

我们的第二个方法是极端随机化树（Extremely Randomized trees），这是一种与随机森林非常相似的算法。通过与随机森林相比，增加更多的随机化，它声称可以更有效地解决方差问题。此外，它还稍微减少了计算复杂度。

我们的最后一个方法是旋转森林（Rotation Forest）。前两个方法需要大量的树作为其集成的一部分，以获得良好的性能。旋转森林声称可以用较少的树实现类似或更好的性能。此外，该算法的作者声称，其基础估计器可以是任何其他的模型，而不仅仅是树。通过这种方式，它被视为构建类似于梯度提升（Gradient Boosting）集成的新框架。

从树到森林——随机森林

随机森林方法构建了许多相互之间不相关的树（森林）。给定一个分类或回归问题，该方法构建许多树，最终的预测结果要么是森林中所有树的预测平均值（对于回归），要么是多数投票分类的结果。

这应该让你想起 Bagging。随机森林是另一种 Bagging 方法。Bagging 背后的基本思想是使用大量的噪声估计器，通过平均来处理噪声，从而减少最终输出中的方差。树对训练数据集中的噪声非常敏感。由于树是噪声估计器，它们非常适合用于 Bagging。

让我们写下构建随机森林的步骤。森林中所需的树的数量是用户指定的一个参数。假设 T 是需要构建的树的数量：

我们从 1 到 T 进行迭代，也就是说，我们构建 T 棵树：

• 对于每棵树，从我们的输入数据集中抽取大小为 D 的自助抽样。

• 我们继续将一棵树 t 拟合到输入数据：

  • 随机选择 m 个属性。

  • 选择最好的属性作为分裂变量，使用预定义的标准。

  • 将数据集分成两部分。记住，树是二叉的。在树的每一层，输入数据集都会被分成两部分。

  • 我们继续在已分割的数据集上递归地执行前面三步。

• 最后，我们返回 T 棵树。

为了对一个新实例做出预测，我们在 T 中所有的树上进行多数投票来做分类；对于回归问题，我们取每棵树 t 在 T 中返回的平均值。

我们之前提到过，随机森林构建的是非相关的树。让我们看看集成中的各个树是如何彼此不相关的。通过为每棵树从数据集中抽取自助样本，我们确保不同的树会接触到数据的不同部分。这样，每棵树都会尝试建模数据集的不同特征。因此，我们遵循集成方法引入底层估计器的变化。但这并不保证底层树之间完全没有相关性。当我们进行节点分裂时，并不是选择所有特征，而是随机选择特征的一个子集。通过这种方式，我们尝试确保我们的树之间没有相关性。

与 Boosting 相比，我们的估计器集成在 Boosting 中是弱分类器，而在随机森林中，我们构建具有最大深度的树，以使其完美拟合自助样本，从而降低偏差。其结果是引入了高方差。然而，通过构建大量的树并使用平均化原则进行最终预测，我们希望解决这个方差问题。

让我们继续深入了解我们的随机森林配方。

准备开始

我们将生成一些分类数据集，以演示随机森林算法。我们将利用 scikit-learn 中的随机森林实现，该实现来自集成模块。

如何操作…

我们将从加载所有必要的库开始。让我们利用 sklearn.dataset 模块中的 make_classification 方法来生成训练数据，以演示随机森林：

from sklearn.datasets import make_classification
from sklearn.metrics import classification_report, accuracy_score
from sklearn.cross_validation import train_test_split
from sklearn.ensemble import RandomForestClassifier
from sklearn.grid_search import RandomizedSearchCV
from operator import itemgetter

import numpy as np

def get_data():
    """
    Make a sample classification dataset
    Returns : Independent variable y, dependent variable x
    """
    no_features = 30
    redundant_features = int(0.1*no_features)
    informative_features = int(0.6*no_features)
    repeated_features = int(0.1*no_features)
    x,y = make_classification(n_samples=500,n_features=no_features,flip_y=0.03,\
            n_informative = informative_features, n_redundant = redundant_features \
            ,n_repeated = repeated_features,random_state=7)
    return x,y

现在我们将编写 build_forest 函数来构建完全生长的树，并继续评估森林的性能。接着我们将编写可用于搜索森林最优参数的方法：

def build_forest(x,y,x_dev,y_dev):
    """
    Build a random forest of fully grown trees
    and evaluate peformance
    """
    no_trees = 100
    estimator = RandomForestClassifier(n_estimators=no_trees)
    estimator.fit(x,y)

    train_predcited = estimator.predict(x)
    train_score = accuracy_score(y,train_predcited)
    dev_predicted = estimator.predict(x_dev)
    dev_score = accuracy_score(y_dev,dev_predicted)

    print "Training Accuracy = %0.2f Dev Accuracy = %0.2f"%(train_score,dev_score)

def search_parameters(x,y,x_dev,y_dev):
    """
    Search the parameters of random forest algorithm
    """
    estimator = RandomForestClassifier()
    no_features = x.shape[1]
    no_iterations = 20
    sqr_no_features = int(np.sqrt(no_features))

    parameters = {"n_estimators"      : np.random.randint(75,200,no_iterations),
                 "criterion"         : ["gini", "entropy"],
                 "max_features"      : [sqr_no_features,sqr_no_features*2,sqr_no_features*3,sqr_no_features+10]
                 }

    grid = RandomizedSearchCV(estimator=estimator,param_distributions=parameters,\
    verbose=1, n_iter=no_iterations,random_state=77,n_jobs=-1,cv=5)
    grid.fit(x,y)
    print_model_worth(grid,x_dev,y_dev)

    return grid.best_estimator_

def print_model_worth(grid,x_dev,y_dev):    
    # Print the goodness of the models
    # We take the top 5 models
    scores = sorted(grid.grid_scores_, key=itemgetter(1), reverse=True) [0:5]

    for model_no,score in enumerate(scores):
        print "Model %d, Score = %0.3f"%(model_no+1,score.mean_validation_score)
        print "Parameters = {0}".format(score.parameters)
    print
    dev_predicted = grid.predict(x_dev)

    print classification_report(y_dev,dev_predicted)

最后，我们编写一个主函数，用于调用我们之前定义的函数：

if __name__ == "__main__":
    x,y = get_data()    

    # Divide the data into Train, dev and test    
    x_train,x_test_all,y_train,y_test_all = train_test_split(x,y,test_size = 0.3,random_state=9)
    x_dev,x_test,y_dev,y_test = train_test_split(x_test_all,y_test_all,test_size=0.3,random_state=9)

    build_forest(x_train,y_train,x_dev,y_dev)
    model = search_parameters(x,y,x_dev,y_dev)
    get_feature_importance(model)

它是如何工作的…

让我们从我们的主函数开始。我们调用 get_data 来获取预测器特征 x 和响应特征 y。在 get_data 中，我们利用 make_classification 数据集来生成我们的随机森林训练数据：

def get_data():
    """
    Make a sample classification dataset
    Returns : Independent variable y, dependent variable x
    """
    no_features = 30
    redundant_features = int(0.1*no_features)
    informative_features = int(0.6*no_features)
    repeated_features = int(0.1*no_features)
    x,y = make_classification(n_samples=500,n_features=no_features,flip_y=0.03,\
            n_informative = informative_features, n_redundant = redundant_features \
            ,n_repeated = repeated_features,random_state=7)
    return x,y

让我们看看传递给make_classification方法的参数。第一个参数是所需的实例数量；在此例中，我们需要 500 个实例。第二个参数是每个实例所需的属性数量。我们设定需要 30 个属性。第三个参数flip_y会随机交换 3%的实例。这是为了向数据中引入一些噪声。接下来的参数指定从这 30 个属性中有多少个是足够有用的信息属性，用于我们的分类任务。我们设定 60%的特征，即 30 个特征中的 18 个应该具有信息量。下一个参数与冗余特征有关。这些冗余特征是通过信息特征的线性组合生成的，以引入特征之间的相关性。最后，重复特征是重复的特征，它们是从信息特征和冗余特征中随机抽取的。

让我们使用train_test_split将数据划分为训练集和测试集。我们将 30%的数据保留用于测试：

    # Divide the data into Train, dev and test    
    x_train,x_test_all,y_train,y_test_all = train_test_split(x,y,test_size = 0.3,random_state=9)

我们再次使用train_test_split将我们的测试数据分成开发集和测试集：

    x_dev,x_test,y_dev,y_test = train_test_split(x_test_all,y_test_all,test_size=0.3,random_state=9)

在数据被分配用于构建、评估和测试模型后，我们继续构建我们的模型：

build_forest(x_train,y_train,x_dev,y_dev)

我们使用训练集和开发集数据调用build_forest函数来构建随机森林模型。让我们来看一下这个函数的内部实现：

    no_trees = 100
    estimator = RandomForestClassifier(n_estimators=no_trees)
    estimator.fit(x,y)

    train_predcited = estimator.predict(x)
    train_score = accuracy_score(y,train_predcited)
    dev_predicted = estimator.predict(x_dev)
    dev_score = accuracy_score(y_dev,dev_predicted)

    print "Training Accuracy = %0.2f Dev Accuracy = %0.2f"%(train_score,dev_score)

我们的集成中需要 100 棵树，所以我们使用变量no_trees来定义树的数量。我们利用 scikit-learn 中的RandomForestClassifier类进行检查并应用。正如你所见，我们将所需的树的数量作为参数传递。然后，我们继续拟合我们的模型。

现在让我们找到我们训练集和开发集的模型准确度分数：

不错！我们在开发集上达到了 83%的准确率。让我们看看是否可以提高我们的得分。随机森林中还有其他可调的参数，可以调节以获得更好的模型。有关可以调节的参数列表，请参考以下链接：

scikit-learn.org/stable/modules/generated/sklearn.ensemble.RandomForestClassifier.html

我们调用search_parameters函数，并使用训练数据和开发数据来调整我们随机森林模型的各项参数。

在一些前面的实例中，我们使用了 GridSearchCV 来遍历参数空间，以寻找最佳的参数组合。GridSearchCV 进行的是非常彻底的搜索。然而，在本实例中，我们将使用 RandomizedSearchCV。我们为每个参数提供一个参数值的分布，并指定所需的迭代次数。在每次迭代中，RandomizedSearchCV 将从参数分布中随机选择一个值并拟合模型：

parameters = {"n_estimators" : np.random.randint(75,200,no_iterations),
"criterion" : ["gini", "entropy"],
"max_features" : [sqr_no_features,sqr_no_features*2,sqr_no_features*3,sqr_no_features+10]
}

我们提供一个参数字典，就像我们在 GridSearchCV 中做的那样。在我们的情况下，我们想要测试三个参数。

第一个参数是模型中的树木数量，通过n_estimators参数表示。通过调用 randint 函数，我们获得一个 75 到 200 之间的整数列表。树木的大小由no_iterations参数定义：

no_iterations = 20

这是我们将传递给 RandomizedSearchCV 的参数，表示我们希望执行的迭代次数。在这20个元素的数组中，RandomizedSearchCV 将为每次迭代随机抽取一个值。

下一个参数是准则，我们在基尼指数和熵之间随机选择，并将其作为每次迭代中拆分节点的准则。

最重要的参数max_features定义了算法在拆分每个节点时应该选择的特征数量。在我们描述随机森林的伪代码中，我们指定了每次拆分节点时需要随机选择 m 个特征。max_features参数定义了 m。在这里，我们提供了一个包含四个值的列表。变量sqr_no_features是输入数据集中特征数量的平方根：

sqr_no_features = int(np.sqrt(no_features))

列表中的其他值是平方根的一些变化。

让我们用这个参数分布来实例化 RandomizedSearchCV：

grid = RandomizedSearchCV(estimator=estimator,param_distributions=parameters,\
verbose=1, n_iter=no_iterations,random_state=77,n_jobs=-1,cv=5)

第一个参数是底层估算器，即我们试图优化其参数的模型。它是我们的RandomForestClassifier：

estimator = RandomForestClassifier()

第二个参数param_distributions是通过字典参数定义的分布。我们定义了迭代次数，即我们希望运行 RandomForestClassifier 的次数，使用参数n_iter。通过cv参数，我们指定所需的交叉验证次数，在我们的例子中是5次交叉验证。

让我们继续拟合模型，看看模型的效果如何：

grid.fit(x,y)
print_model_worth(grid,x_dev,y_dev)

如你所见，我们有五个折叠，也就是说，我们希望在每次迭代中进行五折交叉验证。我们总共执行20次迭代，因此我们将构建 100 个模型。

让我们看看print_model_worth函数内部。我们将网格对象和开发数据集传递给这个函数。网格对象在一个名为grid_scores_的属性中存储了它构建的每个模型的评估指标，这个属性是一个列表。让我们将这个列表按降序排序，以构建最佳模型：

scores = sorted(grid.grid_scores_, key=itemgetter(1), reverse=True) [0:5]

我们选择排名前五的模型，如索引所示。接下来我们打印这些模型的详细信息：

for model_no,score in enumerate(scores):
print "Model %d, Score = %0.3f"%(model_no+1,score.mean_validation_score)
print "Parameters = {0}".format(score.parameters)
print

我们首先打印评估分数，并接着展示模型的参数：

我们根据模型的得分将模式按降序排列，从而将最佳的模型参数放在最前面。我们将选择这些参数作为我们的模型参数。属性best_estimator_将返回具有这些参数的模型。

让我们使用这些参数并测试我们的开发数据：

dev_predicted = grid.predict(x_dev)
print classification_report(y_dev,dev_predicted)

predict 函数将内部使用best_estimtor：

太棒了！我们有一个完美的模型，分类准确率为 100%。

还有更多…

在内部，RandomForestClassifier使用DecisionTreeClassifier。有关构建决策树时传递的所有参数，请参考以下链接：

scikit-learn.org/stable/modules/generated/sklearn.tree.DecisionTreeClassifier.html

一个我们感兴趣的参数是 splitter。splitter 的默认值设置为 best。根据max_features属性，内部实现将选择划分机制。可用的划分机制包括以下几种：

• best: 从由max_features参数定义的给定属性集中选择最佳可能的划分

• random: 随机选择一个划分属性

你可能已经注意到，在实例化RandomForestClassifier时，这个参数不可用。唯一的控制方式是为max_features参数赋值，该值应小于数据集中的属性数量。

在工业界，随机森林被广泛用于变量选择。在 Scikit learn 中，变量重要性是通过基尼不纯度（gini impurity）来计算的。用于节点划分的基尼和熵准则根据它们将数据集划分为高不纯度的子集的能力来识别最佳划分属性，以便后续的划分能产生良好的分类。一个变量的重要性由它能在划分后的数据集中引入的不纯度量决定。有关更多详细信息，请参考以下书籍：

Breiman, Friedman, "分类与回归树"，1984 年。

我们可以编写一个小函数来打印重要特征：

def get_feature_importance(model):
    feature_importance = model.feature_importances_
    fm_with_id = [(i,importance) for i,importance in enumerate(feature_importance)]
    fm_with_id = sorted(fm_with_id, key=itemgetter(1),reverse=True)[0:10]
    print "Top 10 Features"
    for importance in fm_with_id:
        print "Feature %d importance = %0.3f"%(importance[0],importance[1])
    print

一个 Random Forest 对象有一个名为feature_importances_的变量。我们使用这个变量并创建一个包含特征编号和重要性的元组列表：

    feature_importance = model.feature_importances_
    fm_with_id = [(i,importance) for i,importance in enumerate(feature_importance)]

我们接着按重要性降序排列，并选择前 10 个特征：

    fm_with_id = sorted(fm_with_id, key=itemgetter(1),reverse=True)[0:10]

然后我们打印出前 10 个特征：

随机森林的另一个有趣方面是袋外估计（OOB）。记住，我们最初从数据集中对每棵树进行自助采样（bootstrap）。由于自助采样，某些记录在某些树中不会被使用。假设记录 1 在 100 棵树中被使用，在 150 棵树中未被使用。然后，我们可以使用那 150 棵树来预测该记录的类别标签，从而计算该记录的分类误差。袋外估计可以有效地评估我们森林的质量。以下网址给出了 OOB 如何有效使用的示例：

scikit-learn.org/dev/auto_examples/ensemble/plot_ensemble_oob.html

Scikit learn 中的 RandomForestClassifier 类来源于ForestClassifier。其源代码可以在以下链接找到：

github.com/scikit-learn/scikit-learn/blob/a95203b/sklearn/ensemble/forest.py#L318

当我们在 RandomForestClassifier 中调用predict方法时，它内部调用了在 ForestClassifier 中定义的predict_proba方法。在这里，最终的预测不是通过投票完成，而是通过对森林中不同树的每个类别的概率进行平均，并基于最高概率决定最终类别。

Leo Breiman 关于随机森林的原始论文可以在以下链接下载：

link.springer.com/article/10.1023%2FA%3A1010933404324

你还可以参考 Leo Breiman 和 Adele Cutler 维护的网站：

www.stat.berkeley.edu/~breiman/RandomForests/cc_home.htm

另请参见

• 第六章中的构建决策树以解决多类问题配方，机器学习 I。

• 第八章中的理解集成，梯度提升配方，模型选择与评估。

• 第八章中的理解集成，袋装方法配方，模型选择与评估。

生长极端随机树

极端随机树，也称为 Extra Trees 算法，与前面配方中描述的随机森林在两方面有所不同：

1. 它不使用自助法（bootstrapping）为集成中的每棵树选择实例；相反，它使用完整的训练数据集。

2. 给定 K 作为在某个节点上要随机选择的属性数量，它选择一个随机切分点，而不考虑目标变量。

如前面的配方所示，随机森林在两个地方使用了随机化。首先，选择用于训练森林中树的实例时使用了自助法来选择训练实例。其次，在每个节点，随机选择了一组属性，从中选出一个属性，依据的是基尼不纯度或熵准则。极端随机树更进一步，随机选择分裂属性。

极端随机树在以下论文中被提出：

P. Geurts, D. Ernst., 和 L. Wehenkel，“极端随机树”，《机器学习》，63(1)，3-42，2006。

根据本文，除了之前列出的技术方面，还有两个方面使得极端随机树更为适用：

Extra-Trees 方法背后的原理是，通过显式地随机化切分点和属性，结合集成平均法，应该能够比其他方法使用的较弱随机化方案更强烈地减少方差。

与随机森林相比，切分点的随机化（即在每个节点选择用于切分数据集的属性）结合切分点的随机化，即忽略任何标准，最后平均每棵树的结果，将在未知数据集上表现出更优的性能。

第二个优点是计算复杂度：

从计算角度来看，假设树是平衡的，树的生长过程的复杂度与学习样本大小呈 N log N 的数量级，就像大多数树生长过程一样。然而，考虑到节点切分过程的简洁性，我们预计常数因子将比其他集成方法中局部优化切分点的情况要小得多。

由于没有计算时间用于识别最佳的切分属性，这种方法比随机森林在计算上更高效。

让我们写下构建极度随机化树的步骤。森林中所需的树的数量通常由用户指定。设 T 为需要构建的树的数量。

我们从 1 迭代到 T，也就是我们构建 T 棵树：

• 对于每棵树，我们选择完整的输入数据集。

• 然后我们继续拟合一棵树 t 到输入数据：

  • 随机选择 m 个属性。

  • 随机选择一个属性作为切分变量。

  • 将数据集分成两部分。请记住，树是二叉的。在树的每一层，输入数据集都会被分成两部分。

  • 对我们分割的数据集递归执行前述三个步骤。

• 最后，我们返回 T 棵树。

• 让我们看看极度随机化树的配方。

准备好了...

我们将生成一些分类数据集来演示极度随机化树。为此，我们将利用 Scikit Learn 中极度随机化树集成模块的实现。

如何做到这一点...

我们从加载所有必要的库开始。让我们利用sklearn.dataset模块中的make_classification方法来生成训练数据：

from sklearn.datasets import make_classification
from sklearn.metrics import classification_report, accuracy_score
from sklearn.cross_validation import train_test_split, cross_val_score
from sklearn.ensemble import ExtraTreesClassifier
from sklearn.grid_search import RandomizedSearchCV
from operator import itemgetter

def get_data():
 """
 Make a sample classification dataset
 Returns : Independent variable y, dependent variable x
 """
 no_features = 30
 redundant_features = int(0.1*no_features)
 informative_features = int(0.6*no_features)
 repeated_features = int(0.1*no_features)
 x,y = make_classification(n_samples=500,n_features=no_features,flip_y=0.03,\
 n_informative = informative_features, n_redundant = redundant_features \
 ,n_repeated = repeated_features,random_state=7)
    return x,y

我们编写build_forest函数，在其中构建完全生长的树，并继续评估森林的性能：

def build_forest(x,y,x_dev,y_dev):
    """
    Build a Extremely random tress
    and evaluate peformance
    """
    no_trees = 100
    estimator = ExtraTreesClassifier(n_estimators=no_trees,random_state=51)
    estimator.fit(x,y)

    train_predcited = estimator.predict(x)
    train_score = accuracy_score(y,train_predcited)
    dev_predicted = estimator.predict(x_dev)
    dev_score = accuracy_score(y_dev,dev_predicted)

    print "Training Accuracy = %0.2f Dev Accuracy = %0.2f"%(train_score,dev_score)
    print "cross validated score"
    print cross_val_score(estimator,x_dev,y_dev,cv=5)

def search_parameters(x,y,x_dev,y_dev):
    """
    Search the parameters 
    """
    estimator = ExtraTreesClassifier()
    no_features = x.shape[1]
    no_iterations = 20
    sqr_no_features = int(np.sqrt(no_features))

    parameters = {"n_estimators"      : np.random.randint(75,200,no_iterations),
                 "criterion"         : ["gini", "entropy"],
                 "max_features"      : [sqr_no_features,sqr_no_features*2,sqr_no_features*3,sqr_no_features+10]
                 }

    grid = RandomizedSearchCV(estimator=estimator,param_distributions=parameters,\
    verbose=1, n_iter=no_iterations,random_state=77,n_jobs=-1,cv=5)
    grid.fit(x,y)
    print_model_worth(grid,x_dev,y_dev)

    return grid.best_estimator_

最后，我们编写一个主函数来调用我们定义的函数：

if __name__ == "__main__":
    x,y = get_data()    

    # Divide the data into Train, dev and test    
    x_train,x_test_all,y_train,y_test_all = train_test_split(x,y,test_size = 0.3,random_state=9)
    x_dev,x_test,y_dev,y_test = train_test_split(x_test_all,y_test_all,test_size=0.3,random_state=9)

    build_forest(x_train,y_train,x_dev,y_dev)
    model = search_parameters(x,y,x_dev,y_dev)

它是如何工作的…

让我们从主函数开始。我们调用get_data来获取预测属性和响应属性。在get_data函数内部，我们利用 make_classification 数据集来生成我们配方的训练数据，具体如下：

def get_data():
    """
    Make a sample classification dataset
    Returns : Independent variable y, dependent variable x
    """
    no_features = 30
    redundant_features = int(0.1*no_features)
    informative_features = int(0.6*no_features)
    repeated_features = int(0.1*no_features)
    x,y = make_classification(n_samples=500,n_features=no_features,flip_y=0.03,\
            n_informative = informative_features, n_redundant = redundant_features \
            ,n_repeated = repeated_features,random_state=7)
    return x,y

让我们看看传递给make_classification方法的参数。第一个参数是所需的实例数量；在这个例子中，我们需要 500 个实例。第二个参数是每个实例所需的属性数量。我们需要 30 个属性。第三个参数，flip_y，会随机交换 3%的实例。这样做是为了给我们的数据引入一些噪声。下一个参数指定了这些 30 个特征中，有多少个特征应足够信息量，以便用于分类。我们指定 60%的特征，也就是 30 个中的 18 个，应该是有信息量的。下一个参数是冗余特征。这些特征是通过有信息量的特征的线性组合生成的，以便在特征之间引入相关性。最后，重复特征是从有信息量特征和冗余特征中随机选取的重复特征。

让我们使用train_test_split将数据分为训练集和测试集。我们将 30%的数据用于测试：

    # Divide the data into Train, dev and test    
    x_train,x_test_all,y_train,y_test_all = train_test_split(x,y,test_size = 0.3,random_state=9)

再次，我们使用train_test_split将我们的测试数据分为开发集和测试集：

    x_dev,x_test,y_dev,y_test = train_test_split(x_test_all,y_test_all,test_size=0.3,random_state=9)

数据已经分割用于构建、评估和测试模型，我们继续构建我们的模型：

build_forest(x_train,y_train,x_dev,y_dev)

我们调用build_forest函数，使用我们的训练集和开发集数据来构建极端随机化树模型。让我们来看一下这个函数：

no_trees = 100
    estimator = ExtraTreesClassifier(n_estimators=no_trees,random_state=51)
    estimator.fit(x,y)

    train_predcited = estimator.predict(x)
    train_score = accuracy_score(y,train_predcited)
    dev_predicted = estimator.predict(x_dev)
    dev_score = accuracy_score(y_dev,dev_predicted)

    print "Training Accuracy = %0.2f Dev Accuracy = %0.2f"%(train_score,dev_score)
    print "cross validated score"
    print cross_val_score(estimator,x_dev,y_dev,cv=5)

我们需要 100 棵树来构建我们的集成模型，因此我们使用变量 no_trees 来定义树的数量。我们利用 Scikit Learn 中的ExtraTreesClassifier类。如你所见，我们将所需的树的数量作为参数传递。这里有一个需要注意的点是参数 bootstrap。有关 ExtraTreesClassifier 的参数，请参考以下网址：

scikit-learn.org/stable/modules/generated/sklearn.ensemble.ExtraTreesClassifier.html

参数 bootstrap 默认设置为False。与以下网址给出的RandomForestClassifier的 bootstrap 参数进行对比：

scikit-learn.org/stable/modules/generated/sklearn.ensemble.RandomForestClassifier.html

如前所述，森林中的每棵树都是用所有记录进行训练的。

我们继续如下拟合我们的模型：

    train_predcited = estimator.predict(x)

然后我们继续找到训练集和开发集数据的模型准确度得分：

    train_score = accuracy_score(y,train_predcited)
    dev_predicted = estimator.predict(x_dev)
    dev_score = accuracy_score(y_dev,dev_predicted)

让我们打印出训练集和开发集的数据得分：

    print "Training Accuracy = %0.2f Dev Accuracy = %0.2f"%(train_score,dev_score)

现在让我们进行五折交叉验证，查看模型预测的结果：

相当不错的结果。我们几乎在其中一个折叠中达到了 90%的准确率。我们可以像对随机森林一样在参数空间内进行随机搜索。让我们调用 search_parameters 函数，传入我们的训练和测试数据集。请参阅前面的内容了解 RandomizedSearchCV 的解释。然后，我们将打印 search_parameters 函数的输出：

如同前面的步骤，我们已按得分从高到低对模型进行排名，因此最好的模型参数会排在最前面。我们将选择这些参数作为我们的模型参数。属性 best_estimator_ 将返回具有这些参数的模型。

接下来，您将看到为最佳估计器生成的分类报告。预测函数将内部使用 best_estimator_。报告是通过以下代码生成的：

dev_predicted = grid.predict(x_dev)
print classification_report(y_dev,dev_predicted)

很棒！我们有一个完美的模型，分类准确率为 100%。

还有更多内容……

极端随机树在时间序列分类问题中非常受欢迎。有关更多信息，请参阅以下论文：

Geurts, P., Blanco Cuesta A., 和 Wehenkel, L. (2005a). 生物序列分类的分段与组合方法。收录于：IEEE 生物信息学与计算生物学计算智能研讨会论文集，194–201。

参见

• 构建决策树解决多类别问题 章节内容参见 第六章，机器学习 I

• 理解集成方法，袋装法 章节内容参见 第八章，模型选择与评估

• 从树到森林的生长，随机森林 章节内容参见 第九章，机器学习 III

旋转森林的生长

随机森林和袋装法在非常大的集成中能够给出令人印象深刻的结果；拥有大量估计器会提高这些方法的准确性。相反，旋转森林是设计用来处理较小数量的集成。

让我们写下构建旋转森林的步骤。森林中所需的树木数量通常由用户指定。让 T 表示需要构建的树的数量。

我们从 1 开始迭代直到 T，也就是说，我们构建 T 棵树。

对于每棵树 t，执行以下步骤：

• 将训练集中的特征分成 K 个大小相等的不重叠子集。

• 我们有 K 个数据集，每个数据集包含 K 个特征。对于每个 K 数据集，我们进行如下操作：从每个 K 数据集中自助抽取 75%的数据，并使用抽取的样本进行后续步骤：

  • 对 K 中的第 i 个子集进行主成分分析。保留所有主成分。对于 Kth 子集中的每个特征 j，我们有一个主成分 a。我们将其表示为 aij，这是第 i 个子集中的第 j 个属性的主成分。

  • 存储子集的主成分。

• 创建一个大小为 n X n 的旋转矩阵，其中 n 是属性的总数。将主成分排列在矩阵中，使得这些成分匹配原始训练数据集中特征的位置。

• 使用矩阵乘法将训练数据集投影到旋转矩阵上。

• 使用投影数据集构建决策树。

• 存储树和旋转矩阵。

有了这些知识，让我们跳入我们的配方。

准备工作…

我们将生成一些分类数据集来演示旋转森林。根据我们的了解，目前没有现成的 Python 实现来支持旋转森林。因此，我们将编写自己的代码。我们将利用 Scikit Learn 实现的决策树分类器，并使用train_test_split方法进行自助抽样。

如何实现…

我们将从加载所有必要的库开始。让我们利用来自sklearn.dataset模块的make_classification方法生成训练数据。接着我们使用一个方法来选择一个随机属性子集，称为gen_random_subset：

from sklearn.datasets import make_classification
from sklearn.metrics import classification_report
from sklearn.cross_validation import train_test_split
from sklearn.decomposition import PCA
from sklearn.tree import DecisionTreeClassifier
import numpy as np

def get_data():
    """
    Make a sample classification dataset
    Returns : Independent variable y, dependent variable x
    """
    no_features = 50
    redundant_features = int(0.1*no_features)
    informative_features = int(0.6*no_features)
    repeated_features = int(0.1*no_features)
    x,y = make_classification(n_samples=500,n_features=no_features,flip_y=0.03,\
            n_informative = informative_features, n_redundant = redundant_features \
            ,n_repeated = repeated_features,random_state=7)
    return x,y

def get_random_subset(iterable,k):
    subsets = []
    iteration = 0
    np.random.shuffle(iterable)
    subset = 0
    limit = len(iterable)/k
    while iteration < limit:
        if k <= len(iterable):
            subset = k
        else:
            subset = len(iterable)
        subsets.append(iterable[-subset:])
        del iterable[-subset:]
        iteration+=1
    return subsets

现在我们编写一个函数build_rotationtree_model，在其中我们将构建完全生长的树，并使用model_worth函数评估森林的表现：

def build_rotationtree_model(x_train,y_train,d,k):
    models = []
    r_matrices = []
    feature_subsets = []
    for i in range(d):
        x,_,_,_ = train_test_split(x_train,y_train,test_size=0.3,random_state=7)
        # Features ids
        feature_index = range(x.shape[1])
        # Get subsets of features
        random_k_subset = get_random_subset(feature_index,k)
        feature_subsets.append(random_k_subset)
        # Rotation matrix
        R_matrix = np.zeros((x.shape[1],x.shape[1]),dtype=float)
        for each_subset in random_k_subset:
            pca = PCA()
            x_subset = x[:,each_subset]
            pca.fit(x_subset)
            for ii in range(0,len(pca.components_)):
                for jj in range(0,len(pca.components_)):
                    R_matrix[each_subset[ii],each_subset[jj]] = pca.components_[ii,jj]

        x_transformed = x_train.dot(R_matrix)

        model = DecisionTreeClassifier()
        model.fit(x_transformed,y_train)
        models.append(model)
        r_matrices.append(R_matrix)
    return models,r_matrices,feature_subsets

def model_worth(models,r_matrices,x,y):

    predicted_ys = []
    for i,model in enumerate(models):
        x_mod =  x.dot(r_matrices[i])
        predicted_y = model.predict(x_mod)
        predicted_ys.append(predicted_y)

    predicted_matrix = np.asmatrix(predicted_ys)
    final_prediction = []
    for i in range(len(y)):
        pred_from_all_models = np.ravel(predicted_matrix[:,i])
        non_zero_pred = np.nonzero(pred_from_all_models)[0]  
        is_one = len(non_zero_pred) > len(models)/2
        final_prediction.append(is_one)

    print classification_report(y, final_prediction)

最后，我们编写一个主函数，用于调用我们之前定义的函数：

if __name__ == "__main__":
    x,y = get_data()    
#    plot_data(x,y)

    # Divide the data into Train, dev and test    
    x_train,x_test_all,y_train,y_test_all = train_test_split(x,y,test_size = 0.3,random_state=9)
    x_dev,x_test,y_dev,y_test = train_test_split(x_test_all,y_test_all,test_size=0.3,random_state=9)

    # Build a bag of models
    models,r_matrices,features = build_rotationtree_model(x_train,y_train,25,5)
    model_worth(models,r_matrices,x_train,y_train)
    model_worth(models,r_matrices,x_dev,y_dev)

它是如何工作的…

让我们从主函数开始。我们调用get_data来获取响应属性中的预测器属性。在get_data内部，我们利用make_classification数据集生成我们的训练数据，具体如下：

def get_data():
    """
    Make a sample classification dataset
    Returns : Independent variable y, dependent variable x
    """
    no_features = 30
    redundant_features = int(0.1*no_features)
    informative_features = int(0.6*no_features)
    repeated_features = int(0.1*no_features)
    x,y = make_classification(n_samples=500,n_features=no_features,flip_y=0.03,\
            n_informative = informative_features, n_redundant = redundant_features \
            ,n_repeated = repeated_features,random_state=7)
    return x,y

让我们看一下传递给make_classification方法的参数。第一个参数是所需实例的数量；在这种情况下，我们需要 500 个实例。第二个参数是每个实例所需的属性数量。我们需要 30 个属性。第三个参数flip_y随机交换 3%的实例。这样做是为了在数据中引入一些噪声。下一个参数是关于 30 个特征中应该具有足够信息量用于分类的特征数量。我们规定，60%的特征，也就是 30 个中的 18 个应该具有信息量。下一个参数是冗余特征。冗余特征是通过信息性特征的线性组合生成的，用于在特征之间引入相关性。最后，重复特征是从信息性特征和冗余特征中随机选择的重复特征。

让我们使用train_test_split将数据分割成训练集和测试集。我们将 30%的数据用于测试：

    # Divide the data into Train, dev and test    
    x_train,x_test_all,y_train,y_test_all = train_test_split(x,y,test_size = 0.3,random_state=9)

我们再次利用train_test_split将测试数据分成开发集和测试集，具体如下：

    x_dev,x_test,y_dev,y_test = train_test_split(x_test_all,y_test_all,test_size=0.3,random_state=9)

在数据被划分为用于构建、评估和测试模型的部分后，我们开始构建我们的模型：

    models,r_matrices,features = build_rotationtree_model(x_train,y_train,25,5)

我们调用build_rotationtree_model函数来构建我们的旋转森林。我们传入我们的训练数据、预测变量x_train和响应变量y_train，要构建的树的总数（在本例中为25），以及要使用的特征子集（在本例中为5）。

让我们跳到该函数：

    models = []
    r_matrices = []
    feature_subsets = []

我们首先声明三个列表，用于存储每棵决策树、该树的旋转矩阵，以及在该迭代中使用的特征子集。然后我们继续构建我们的集成中的每棵树。

作为第一项工作，我们进行自助法以保留数据的 75％：

        x,_,_,_ = train_test_split(x_train,y_train,test_size=0.3,random_state=7)

我们利用 Scikit learn 中的train_test_split函数进行自助法（bootstrapping）。然后我们按以下方式决定特征子集：

        # Features ids
        feature_index = range(x.shape[1])
        # Get subsets of features
        random_k_subset = get_random_subset(feature_index,k)
        feature_subsets.append(random_k_subset)

get_random_subset函数接受特征索引和所需 k 个子集的数量作为参数，并返回 K 个子集。

在该函数内部，我们会打乱特征索引。特征索引是一个数字数组，从 0 开始，直到训练集中的特征数量：

    np.random.shuffle(iterable)

假设我们有 10 个特征，我们的 k 值是 5，表示我们需要具有 5 个不重叠特征索引的子集；我们需要进行两次迭代。我们将所需的迭代次数存储在 limit 变量中：

    limit = len(iterable)/k
    while iteration < limit:
        if k <= len(iterable):
            subset = k
        else:
            subset = len(iterable)
        iteration+=1

如果我们的所需子集少于总属性数，我们可以继续使用我们可迭代对象中的前 k 个条目。由于我们已经打乱了可迭代对象，我们将在不同的时间返回不同的数量：

        subsets.append(iterable[-subset:])

在选择一个子集后，我们将其从可迭代对象中移除，因为我们需要不重叠的集合：

        del iterable[-subset:]

当所有子集准备好后，我们按以下方式声明我们的旋转矩阵：

        # Rotation matrix
        R_matrix = np.zeros((x.shape[1],x.shape[1]),dtype=float)

如您所见，我们的旋转矩阵的大小是 n x n，其中 n 是我们数据集中的属性数量。您可以看到，我们使用了 shape 属性来声明这个填充了零的矩阵：

        for each_subset in random_k_subset:
            pca = PCA()
            x_subset = x[:,each_subset]
            pca.fit(x_subset)

对于每个只有 K 个特征的数据子集，我们继续进行主成分分析。

我们按如下方式填充我们的旋转矩阵：

            for ii in range(0,len(pca.components_)):
                for jj in range(0,len(pca.components_)):
                    R_matrix[each_subset[ii],each_subset[jj]] = pca.components_[ii,jj]

例如，假设我们的子集中有三个属性，总共有六个属性。为了说明，假设我们的子集是：

2,4,6 and 1,3,5

我们的旋转矩阵 R 是 6 x 6 的大小。假设我们要填充第一个特征子集的旋转矩阵。我们将有三个主成分，分别对应于 2、4 和 6，大小为 1 x 3。

来自 Scikit learn 的 PCA 输出是一个大小为主成分 X 特征的矩阵。我们通过 for 循环遍历每个主成分值。在第一次运行时，我们感兴趣的特征是 2，来自 PCA 的主成分矩阵（0,0）的单元格给出了特征 2 对主成分 1 的贡献值。我们需要找到旋转矩阵中这个值的位置。我们使用主成分矩阵中的索引 ii 和 jj 与子集列表结合来找到旋转矩阵中的正确位置：

                    R_matrix[each_subset[ii],each_subset[jj]] = pca.components_[ii,jj]

each_subset[0]和each_subset[0]将使我们处于旋转矩阵中的(2,2)单元格。随着循环的进行，组件矩阵中(0,1)单元格中的下一个组件值将被放置到旋转矩阵中的(2,4)单元格，最后一个将放在旋转矩阵中的(2,6)单元格。这对于第一个子集中的所有属性都是如此。让我们进入第二个子集；这里第一个属性是 1。组件矩阵中的(0,0)单元格对应于旋转矩阵中的(1,1)单元格。

以这种方式进行，你会发现属性组件的值与属性本身的顺序一致。

在旋转矩阵准备好后，让我们将输入数据投影到旋转矩阵上：

        x_transformed = x_train.dot(R_matrix)

现在是时候拟合我们的决策树了：

model = DecisionTreeClassifier()
model.fit(x_transformed,y_train)

最后，我们存储我们的模型和相应的旋转矩阵：

models.append(model)
r_matrices.append(R_matrix)

构建好模型后，让我们使用model_worth函数来检验模型在训练数据和开发数据上的表现：

model_worth(models,r_matrices,x_train,y_train)
model_worth(models,r_matrices,x_dev,y_dev)

让我们看看我们的model_worth函数：

    for i, model in enumerate(models):
        x_mod =  x.dot(r_matrices[i])
        predicted_y = model.predict(x_mod)
        predicted_ys.append(predicted_y)

在这个函数内部，我们使用每棵树进行预测。然而，在预测之前，我们先使用旋转矩阵对输入数据进行投影。我们将所有预测的输出存储在一个名为predicted_ys的列表中。假设我们有 100 个实例需要预测，并且我们有 10 个模型。在每个实例中，我们会有 10 个预测结果。为了方便起见，我们将这些结果存储为一个矩阵：

predicted_matrix = np.asmatrix(predicted_ys)

现在我们继续为每个输入记录给出最终分类：

    final_prediction = []
    for i in range(len(y)):
        pred_from_all_models = np.ravel(predicted_matrix[:,i])
        non_zero_pred = np.nonzero(pred_from_all_models)[0]  
        is_one = len(non_zero_pred) > len(models)/2
        final_prediction.append(is_one)

我们会将最终的预测结果存储在一个名为final_prediction的列表中。我们遍历每个实例的预测结果。假设我们在第一个实例（在 for 循环中 i=0）；pred_from_all_models存储来自我们模型中所有树的输出。它是一个由 0 和 1 组成的数组，表示在该实例中模型所分类的类别。

我们将其转换为另一个数组non_zero_pred，该数组只包含父数组中非零的条目。

最后，如果这个非零数组的长度大于我们模型数量的一半，我们就说我们最终的预测结果是该实例的类别为 1。我们在这里实现的是经典的投票机制。

让我们通过调用classification_report来看看我们的模型现在有多好：

    print classification_report(y, final_prediction)

以下是我们模型在训练集上的表现：

让我们看看我们的模型在开发数据集上的表现：

还有更多内容…

有关旋转森林的更多信息可以参考以下论文：

旋转森林：一种新的分类器集成方法，Juan J. Rodriguez，IEEE 计算机学会会员，Ludmila I. Kuncheva，IEEE 会员，Carlos J. Alonso

该论文还声称，在 33 个数据集上将极度随机化树与 Bagging、AdBoost 和随机森林进行比较时，极度随机化树的表现优于其他三种算法。

与梯度提升法类似，论文的作者声称极度随机化方法是一个总体框架，且基础集成方法不一定非得是决策树。目前正在进行其他算法的测试工作，如朴素贝叶斯、神经网络等。

另见

• 提取主成分 配方见于第四章，数据分析 - 深度解析

• 通过随机投影减少数据维度 配方见于第四章，数据分析 - 深度解析

• 建立决策树以解决多分类问题 配方见于第六章，机器学习 I

• 理解集成方法，梯度提升法 配方见于第八章，模型选择与评估

• 从树到森林，随机森林 配方见于第九章，机器学习 III

• 生长极度随机化树 配方见于第九章，机器学习 III

第十章 大规模机器学习 – 在线学习

在本章中，我们将看到以下内容：

• 使用感知机作为在线线性算法

• 使用随机梯度下降进行回归

• 使用随机梯度下降进行分类

引言

本章中，我们将专注于大规模机器学习以及适合解决此类大规模问题的算法。直到现在，当我们训练所有模型时，我们假设训练集可以完全存入计算机的内存中。在本章中，我们将看到如何在这一假设不再成立时构建模型。我们的训练数据集非常庞大，因此无法完全加载到内存中。我们可能需要按块加载数据，并且仍能生成一个具有良好准确度的模型。训练集无法完全存入计算机内存的问题可以扩展到流数据上。对于流数据，我们不会一次性看到所有的数据。我们应该能够根据我们所接触到的数据做出决策，并且还应有一个机制，在新数据到达时不断改进我们的模型。

我们将介绍基于随机梯度下降的算法框架。这是一个多功能的框架，用于处理非常大规模的数据集，这些数据集无法完全存入我们的内存中。包括逻辑回归、线性回归和线性支持向量机等多种类型的线性算法，都可以通过该框架进行处理。我们在前一章中介绍的核技巧，也可以纳入该框架中，以应对具有非线性关系的数据集。

我们将从感知机算法开始介绍，感知机是最古老的机器学习算法。感知机易于理解和实现。然而，感知机仅限于解决线性问题。基于核的感知机可以用于解决非线性数据集。

在我们的第二个案例中，我们将正式介绍基于梯度下降的方法框架，以及如何使用该方法执行回归任务。我们将查看不同的损失函数，看看如何使用这些函数构建不同类型的线性模型。我们还将看到感知机如何属于随机梯度下降家族。

在我们的最后一个案例中，我们将看到如何使用随机梯度下降框架构建分类算法。

即使我们没有直接的流数据示例，通过现有的数据集，我们将看到如何解决流数据的使用场景。在线学习算法不限于流数据，它们也可以应用于批量数据，只是它们一次只处理一个实例。

使用感知机作为在线学习算法

如前所述，感知机是最古老的机器学习算法之一。它最早在 1943 年的一篇论文中提到：

《神经活动中固有思想的逻辑演算》。沃伦·S·麦卡洛克和沃尔特·皮茨，伊利诺伊大学医学学院，伊利诺伊神经精神病学研究所精神病学系，芝加哥大学，美国芝加哥。

让我们重新审视分类问题的定义。每个记录或实例可以表示为一个集合(X, y)，其中 X 是一组属性，y 是相应的类标签。

学习一个目标函数 F，该函数将每个记录的属性集映射到预定义的类标签 y，是分类算法的任务。

我们的差异在于我们面对的是一个大规模的学习问题。我们的所有数据无法完全加载到主内存中。因此，我们需要将数据保存在磁盘上，并且每次只使用其中的一部分来构建我们的感知器模型。

让我们继续概述感知器算法：

1. 将模型的权重初始化为一个小的随机数。

2. 将输入数据x以其均值进行中心化。

3. 在每个时间步 t（也称为 epoch）：

   • 洗牌数据集

   • 选择单个记录实例并进行预测。

   • 观察预测与真实标签输出的偏差。

   • 如果预测与真实标签不同，则更新权重。

让我们考虑以下场景。我们将完整的数据集存储在磁盘上。在单个 epoch 中，即在步骤 3 中，所有提到的步骤都在磁盘上的所有数据上执行。在在线学习场景中，一组基于窗口函数的实例将随时提供给我们。我们可以在单个 epoch 中根据窗口中的实例数量更新权重。

让我们来看一下如何更新我们的权重。

假设我们的输入 X 如下：

我们的Y如下：

我们将定义我们的权重如下方程：

我们在看到每个记录后的预测定义如下：

符号函数在权重和属性的乘积为正时返回+1，如果乘积为负，则返回-1。

感知器接着将预测的 y 与实际的 y 进行比较。如果预测的 y 是正确的，它将继续处理下一个记录。如果预测错误，则有两种情况。如果预测的 y 是+1，而实际的 y 是-1，则减小该权重值与 x 的乘积，反之亦然。如果实际的 y 是+1，而预测的 y 是-1，则增加权重。我们可以用方程式来更清楚地表示这一点：

通常，学习率 alpha 会被提供，以便以可控的方式更新权重。由于数据中存在噪声，完全的增量和减量将导致权重无法收敛：

Alpha 是一个非常小的值，范围在 0.1 到 0.4 之间。

现在让我们跳入我们的配方。

准备工作

我们将使用 make_classification 生成数据，采用生成器函数批量生成，以模拟大规模数据和数据流，并继续编写感知机算法。

如何实现…

让我们加载必要的库。然后我们将编写一个名为 get_data 的生成器函数：

from sklearn.datasets import make_classification
from sklearn.metrics import  classification_report
from sklearn.preprocessing import scale
import numpy as np

def get_data(batch_size):
    """
    Make a sample classification dataset
    Returns : Independent variable y, dependent variable x
    """
    b_size = 0
    no_features = 30
    redundant_features = int(0.1*no_features)
    informative_features = int(0.8*no_features)
    repeated_features = int(0.1*no_features)

    while b_size < batch_size:
        x,y = make_classification(n_samples=1000,n_features=no_features,flip_y=0.03,\
                n_informative = informative_features, n_redundant = redundant_features \
                ,n_repeated = repeated_features, random_state=51)
        y_indx = y < 1
        y[y_indx] = -1
        x = scale(x,with_mean=True,with_std=True)

        yield x,y
        b_size+=1

我们将编写两个函数，一个用来构建感知机模型，另一个用来测试我们的模型的有效性：

def build_model(x,y,weights,epochs,alpha=0.5):
    """
    Simple Perceptron
    """

    for i in range(epochs):

        # Shuffle the dataset
        shuff_index = np.random.shuffle(range(len(y)))
        x_train = x[shuff_index,:].reshape(x.shape)
        y_train = np.ravel(y[shuff_index,:])

        # Build weights one instance at a time
        for index in range(len(y)):
            prediction = np.sign( np.sum(x_train[index,:] * weights) ) 
            if prediction != y_train[index]:
                weights = weights + alpha * (y_train[index] * x_train[index,:])

    return weights

def model_worth(x,y,weights):
    prediction = np.sign(np.sum(x * weights,axis=1))
print classification_report(y,prediction)

最后，我们将编写主函数来调用所有前面的函数，以展示感知机算法：

if __name__ == "__main__":
    data = get_data(10)    
    x,y = data.next()
    weights = np.zeros(x.shape[1])    
    for i in range(10):
        epochs = 100
        weights = build_model(x,y,weights,epochs)
        print
        print "Model worth after receiving dataset batch %d"%(i+1)    
        model_worth(x,y,weights)
        print
        if i < 9:
            x,y = data.next()

工作原理…

让我们从主函数开始。我们将要求生成器给我们发送 10 组数据：

    data = get_data(10)    

在这里，我们希望模拟大规模数据和数据流。在构建模型时，我们无法访问所有数据，只能访问部分数据：

    x,y = data.next()

我们将使用生成器中的 next() 函数来获取下一组数据。在 get_data 函数中，我们将使用 scikit-learn 的 make_classification 函数：

        x,y = make_classification(n_samples=1000,n_features=no_features,flip_y=0.03,\
                n_informative = informative_features, n_redundant = redundant_features \
                ,n_repeated = repeated_features, random_state=51)

让我们看看传递给 make_classification 方法的参数。第一个参数是所需的实例数量，在本例中，我们需要 1,000 个实例。第二个参数是每个实例所需的属性数量。我们假设需要 30 个属性。第三个参数 flip_y 随机交换 3% 的实例。这是为了在数据中引入一些噪声。下一个参数涉及这 30 个特征中有多少个应该足够有信息量，以便用于我们的分类。我们指定 60% 的特征，即 30 个特征中的 18 个，应该是有信息量的。下一个参数是冗余特征。这些特征是通过有信息量的特征的线性组合生成的，以便在特征之间引入相关性。最后，重复特征是从有信息量的特征和冗余特征中随机抽取的重复特征。

当我们调用 next() 时，我们将获得 1,000 个数据实例。该函数返回一个 y 标签，值为 {0,1}；我们希望得到 {-1,+1}，因此我们将 y 中的所有零改为 -1：

        y_indx = y < 1
        y[y_indx] = -1

最后，我们将使用 scikit-learn 的 scale 函数来对数据进行中心化处理。

让我们使用第一批数据来构建模型。我们将用零初始化我们的权重矩阵：

    weights = np.zeros(x.shape[1])    

由于我们需要 10 批数据来模拟大规模学习和数据流，因此我们将在 for 循环中执行 10 次模型构建：

    for i in range(10):
        epochs = 100
        weights = build_model(x,y,weights,epochs)

我们的感知机算法在build_model中构建。一个预测变量 x、响应变量 y、权重矩阵和时间步数或周期数作为参数传递。在我们的情况下，我们已将周期数设置为100。此函数还有一个附加参数：alpha 值：

def build_model(x,y,weights,epochs,alpha=0.5)

默认情况下，我们将 alpha 值设置为0.5。

让我们在build_model中查看。我们将从数据洗牌开始：

        # Shuffle the dataset
        shuff_index = np.random.shuffle(range(len(y)))
        x_train = x[shuff_index,:].reshape(x.shape)
        y_train = np.ravel(y[shuff_index,:])

我们将遍历数据集中的每条记录，并开始更新我们的权重：

 # Build weights one instance at a time
        for index in range(len(y)):
            prediction = np.sign( np.sum(x_train[index,:] * weights) ) 
            if prediction != y_train[index]:
                weights = weights + alpha * (y_train[index] * x_train[index,:])

在 for 循环中，你可以看到我们在做预测：

            prediction = np.sign( np.sum(x_train[index,:] * weights) ) 

我们将把训练数据与权重相乘，并将它们加在一起。最后，我们将使用 np.sign 函数来获取我们的预测结果。根据预测结果，我们将更新我们的权重：

                weights = weights + alpha * (y_train[index] * x_train[index,:])

就是这样。我们将把权重返回给调用函数。

在我们的主函数中，我们将调用model_worth函数来打印模型的优度。这里，我们将使用classification_report便捷函数来打印模型的准确度评分：

        print
        print "Model worth after receiving dataset batch %d"%(i+1)    
        model_worth(x,y,weights)

然后我们将继续更新模型以处理下一批输入数据。请注意，我们没有更改weights参数。它会随着每批新数据的到来而更新。

让我们看看model_worth打印了什么：

还有更多…

Scikit-learn 为我们提供了感知机的实现。更多细节请参考以下网址：

scikit-learn.org/stable/modules/generated/sklearn.linear_model.Perceptron.html。

感知机算法中可以改进的另一个方面是使用更多特征。

记住预测方程，我们可以将其重写如下：

我们用一个函数替换了 x 值。在这里，我们可以发送一个特征生成器。例如，一个多项式特征生成器可以添加到我们的get_data函数中，如下所示：

def get_data(batch_size):
    """
    Make a sample classification dataset
    Returns : Independent variable y, dependent variable x
    """
    b_size = 0
    no_features = 30
    redundant_features = int(0.1*no_features)
    informative_features = int(0.8*no_features)
    repeated_features = int(0.1*no_features)
    poly = PolynomialFeatures(degree=2)

    while b_size < batch_size:
        x,y = make_classification(n_samples=1000,n_features=no_features,flip_y=0.03,\
                n_informative = informative_features, n_redundant = redundant_features \
                ,n_repeated = repeated_features, random_state=51)
        y_indx = y < 1
        y[y_indx] = -1
        x = poly.fit_transform(x)
        yield x,y
        b_size+=1

最后，基于核的感知机算法可以处理非线性数据集。有关基于核的感知机的更多信息，请参阅维基百科文章：

en.wikipedia.org/wiki/Kernel_perceptron。

另见

• 在第五章中，学习和使用核的方法，数据挖掘 - 在大海捞针中找针。

使用随机梯度下降进行回归

在典型的回归设置中，我们有一组预测变量（实例），如下所示：

每个实例有 m 个属性，如下所示：

响应变量 Y 是一个实值条目的向量。回归的任务是找到一个函数，使得当 x 作为输入提供给该函数时，它应返回 y：

前面的函数通过一个权重向量进行参数化，也就是说，使用权重向量和输入向量的组合来预测Y，因此通过权重向量重新编写函数会得到如下形式：

所以，现在的问题是我们如何知道我们拥有正确的权重向量？我们将使用损失函数 L 来得到正确的权重向量。损失函数衡量了预测错误的成本。它通过经验来衡量预测 y 时的成本，当实际值为 y 时。回归问题现在变成了寻找正确的权重向量的问题，该向量能够最小化损失函数。对于我们包含n个元素的整个数据集，总的损失函数如下：

我们的权重向量应该是那些能够最小化前面值的权重向量。

梯度下降是一种优化技术，用于最小化前面的方程。对于这个方程，我们将计算梯度，即相对于 W 的一阶导数。

不同于批量梯度下降等其他优化技术，随机梯度下降每次只操作一个实例。随机梯度下降的步骤如下：

1. 对于每一个周期，打乱数据集。

2. 选择一个实例及其响应变量 y。

3. 计算损失函数及其相对于权重的导数。

4. 更新权重。

假设：

这表示相对于 w 的导数。权重更新如下：

如你所见，权重朝着与梯度相反的方向移动，从而实现下降，最终得到的权重向量值能够减少目标成本函数。

平方损失是回归中常用的损失函数。一个实例的平方损失定义如下：

前面方程的导数被代入到权重更新方程中。通过这些背景知识，让我们继续介绍随机梯度下降回归的步骤。

如感知机中所解释的，学习率 eta 被添加到权重更新方程中，以避免噪声的影响：

准备工作

我们将利用 scikit-learn 实现的 SGD 回归。与之前的一些示例一样，我们将使用 scikit-learn 中的make_regression函数生成数据，以展示随机梯度下降回归。

如何操作…

让我们从一个非常简单的示例开始，演示如何构建一个随机梯度下降回归器。

我们首先将加载所需的库。然后我们将编写一个函数来生成预测变量和响应变量，以演示回归：

from sklearn.datasets import make_regression
from sklearn.linear_model import SGDRegressor
from sklearn.metrics import mean_absolute_error, mean_squared_error
from sklearn.cross_validation import train_test_split

def get_data():
    """
    Make a sample classification dataset
    Returns : Independent variable y, dependent variable x
    """
    no_features = 30

    x,y = make_regression(n_samples=1000,n_features=no_features,\
             random_state=51)
    return x,y

我们将继续编写有助于我们构建、验证和检查模型的函数：

def build_model(x,y):
    estimator = SGDRegressor(n_iter = 10, shuffle=True,loss = "squared_loss", \
            learning_rate='constant',eta0=0.01,fit_intercept=True, \
            penalty='none')
    estimator.fit(x,y)

    return estimator

def model_worth(model,x,y):
    predicted_y = model.predict(x)
    print "\nMean absolute error = %0.2f"%mean_absolute_error(y,predicted_y)
    print "Mean squared error = %0.2f"%mean_squared_error(y,predicted_y)

def inspect_model(model):
    print "\nModel Itercept {0}".format(model.intercept_)
    print
    for i,coef in enumerate(model.coef_):
        print "Coefficient {0} = {1:.3f}".format(i+1,coef)

最后，我们将编写我们的主函数来调用所有前面的函数：

if __name__ == "__main__":
    x,y = get_data()

    # Divide the data into Train, dev and test    
    x_train,x_test_all,y_train,y_test_all = train_test_split(x,y,test_size = 0.3,random_state=9)
    x_dev,x_test,y_dev,y_test = train_test_split(x_test_all,y_test_all,test_size=0.3,random_state=9)

    model = build_model(x_train,y_train)

    inspect_model(model)

    print "Model worth on train data"
    model_worth(model,x_train,y_train)
    print "Model worth on dev data"
    model_worth(model,x_dev,y_dev)

    # Building model with l2 regularization
    model = build_model_regularized(x_train,y_train)
    inspect_model(model)

工作原理…

让我们从主函数开始。我们将调用get_data函数来生成我们的预测变量 x 和响应变量 y：

    x,y = get_data()

在get_data函数中，我们将利用 scikit-learn 中便捷的make_regression函数来生成回归问题的数据集：

    no_features = 30
    x,y = make_regression(n_samples=1000,n_features=no_features,\
             random_state=51)	

如您所见，我们将生成一个包含 1,000 个实例的数据集，其中实例数量由n_samples参数指定，特征数量由n_features参数定义，共有 30 个特征。

让我们使用train_test_split将数据划分为训练集和测试集。我们将保留 30%的数据用于测试：

    # Divide the data into Train, dev and test    
    x_train,x_test_all,y_train,y_test_all = train_test_split(x,y,test_size = 0.3,random_state=9)

我们将再次利用train_test_split来将测试数据划分为开发集和测试集：

    x_dev,x_test,y_dev,y_test = train_test_split(x_test_all,y_test_all,test_size=0.3,random_state=9)

在数据划分为用于构建、评估和测试模型之后，我们将继续构建我们的模型。

我们将使用训练数据集调用build_model函数：

    model = build_model(x_train,y_train)

在build_model中，我们将利用 scikit-learn 的 SGD 回归器类来构建我们的随机梯度下降方法：

    estimator = SGDRegressor(n_iter = 10, shuffle=True,loss = "squared_loss", \
            learning_rate='constant',eta0=0.01,fit_intercept=True, \
            penalty='none')
    estimator.fit(x,y)

SGD 回归器是一种广泛的方法，可以用于拟合具有大量参数的多种线性模型。我们将首先解释随机梯度下降的基本方法，然后再继续解释其他细节。

让我们看看我们使用的参数。第一个参数是我们希望通过数据集的次数，用于更新权重。这里，我们将设置为 10 次迭代。与感知机类似，在遍历所有记录一次之后，我们需要在开始下一次迭代时打乱输入记录。我们使用shuffle参数来实现这一点。shuffle的默认值为true，我们这里包含它是为了说明。我们的损失函数是平方损失，并且我们要进行线性回归，因此我们将通过loss参数指定这一点。

我们的学习率eta是一个常数，我们将通过learning_rate参数指定。我们将通过eta=0来为学习率提供一个值。然后我们会说我们需要拟合截距，因为我们没有按数据的均值来对数据进行中心化。最后，惩罚参数控制所需的收缩类型。在我们的案例中，我们不需要任何收缩，因此使用none字符串。

我们将通过调用拟合函数并传入我们的预测变量和响应变量来构建模型。最后，我们将把构建好的模型返回给调用函数。

现在让我们检查一下我们的模型，并查看截距和系数的值：

    inspect_model(model)

在检查模型时，我们将打印模型截距和系数的值：

现在让我们来看一下我们的模型在训练数据中的表现：

    print "Model worth on train data"
    model_worth(model,x_train,y_train)

我们将调用 model_worth 函数来查看模型的表现。model_worth 函数打印出平均绝对误差和均方误差的值。

均方误差定义如下：

平均绝对误差定义如下：

均方误差对异常值非常敏感。因此，平均绝对误差是一个更稳健的衡量标准。让我们通过训练数据来查看模型的表现：

现在让我们来看一下使用开发数据的模型表现：

还有更多…

我们可以在随机梯度下降框架中加入正则化。回顾前一章节中岭回归的成本函数：

在这里，我们加入了平方损失函数的扩展版本，并添加了正则化项——权重平方和。我们可以将其包含在我们的梯度下降过程中。假设我们将正则化项表示为 R(W)。我们的权重更新规则现在如下：

如你所见，现在我们有了损失函数关于权重向量 w 的导数，正则化项对权重的导数被添加到我们的权重更新规则中。

让我们写一个新的函数来构建包含正则化的模型：

def build_model_regularized(x,y):
    estimator = SGDRegressor(n_iter = 10,shuffle=True,loss = "squared_loss", \
            learning_rate='constant',eta0=0.01,fit_intercept=True, \
            penalty='l2',alpha=0.01)
    estimator.fit(x,y)

    return estimator

我们可以通过如下方式从主函数中调用这个函数：

model = build_model_regularized(x_train,y_train)
inspect_model(model)

让我们来看看与之前构建模型方法相比，我们传递的新参数：

    estimator = SGDRegressor(n_iter = 10,shuffle=True,loss = "squared_loss", \
            learning_rate='constant',eta0=0.01,fit_intercept=True, \
            penalty='l2',alpha=0.01)

之前，我们提到过惩罚项为 none。现在，你可以看到我们提到需要在模型中加入 L2 惩罚项。我们将给 alpha 参数一个值为 0.01。让我们来看看我们的系数：

你可以看到 L2 正则化的效果：很多系数已被压缩为零。类似地，L1 正则化和结合了 L1 和 L2 正则化的弹性网方法，也可以通过惩罚参数来实现。

记得在介绍中提到过，随机梯度下降更像是一个框架，而不是单一的方法。通过更改损失函数，可以使用这个框架生成其他线性模型。

可以使用不敏感于 epsilon 的损失函数构建支持向量机回归模型。这个损失函数定义如下：

请参考以下网址，了解可以传递给 scikit-learn 中 SGD 回归器的各种参数：

scikit-learn.org/stable/modules/generated/sklearn.linear_model.SGDRegressor.html。

另见

• 在第七章中，使用回归预测实数值 的方法，机器学习 II

• 第七章中的岭回归的收缩示例，机器学习 II

使用随机梯度下降进行分类

分类问题的设置与回归问题非常相似，唯一的区别在于响应变量。在分类问题中，响应是一个类别变量。由于其性质，我们有不同的损失函数来衡量错误预测的代价。假设我们的讨论和配方是针对二分类器的，我们的目标变量 Y 可以取值 {0,1}。

我们将使用该损失函数的导数作为权重更新规则，以得到我们的权重向量。

scikit-learn 的 SGD 分类器类为我们提供了多种损失函数。然而，在本示例中，我们将看到对数损失，它将给出逻辑回归。

逻辑回归为以下形式的数据拟合一个线性模型：

我们已经给出了一个广义的符号表示。假设截距是我们权重向量的第一维。对于二分类问题，应用对数几率函数来得到预测，如下所示：

上述函数也被称为 sigmoid 函数。对于非常大的正值 x_i，该函数将返回接近 1 的值，反之，对于非常大的负值，将返回接近 0 的值。由此，我们可以将对数损失函数定义如下：

将上述损失函数应用于梯度下降的权重更新规则，我们可以得到适当的权重向量。

对于在 scikit-learn 中定义的对数损失函数，请参考以下网址：

scikit-learn.org/stable/modules/generated/sklearn.metrics.log_loss.html。

了解这些知识后，让我们进入基于随机梯度下降的分类配方。

准备工作

我们将利用 scikit-learn 实现的随机梯度下降分类器。就像在之前的一些示例中一样，我们将使用 scikit-learn 的 make_classification 函数来生成数据，以演示随机梯度下降分类。

如何实现……

我们从一个非常简单的例子开始，演示如何构建一个随机梯度下降回归器。

我们将首先加载所需的库。然后，我们将编写一个函数来生成预测变量和响应变量：

from sklearn.datasets import make_classification
from sklearn.metrics import  accuracy_score
from sklearn.cross_validation import train_test_split
from sklearn.linear_model import SGDClassifier

import numpy as np

def get_data():
    """
    Make a sample classification dataset
    Returns : Independent variable y, dependent variable x
    """
    no_features = 30
    redundant_features = int(0.1*no_features)
    informative_features = int(0.6*no_features)
    repeated_features = int(0.1*no_features)
    x,y = make_classification(n_samples=1000,n_features=no_features,flip_y=0.03,\
            n_informative = informative_features, n_redundant = redundant_features \
            ,n_repeated = repeated_features,random_state=7)
    return x,y

我们将继续编写一些函数，帮助我们构建和验证我们的模型：

def build_model(x,y,x_dev,y_dev):
    estimator = SGDClassifier(n_iter=50,shuffle=True,loss="log", \
                learning_rate = "constant",eta0=0.0001,fit_intercept=True, penalty="none")
    estimator.fit(x,y)
    train_predcited = estimator.predict(x)
    train_score = accuracy_score(y,train_predcited)
    dev_predicted = estimator.predict(x_dev)
    dev_score = accuracy_score(y_dev,dev_predicted)

    print 
    print "Training Accuracy = %0.2f Dev Accuracy = %0.2f"%(train_score,dev_score)

最后，我们将编写主函数，调用所有先前的函数：

if __name__ == "__main__":
    x,y = get_data()    

    # Divide the data into Train, dev and test    
    x_train,x_test_all,y_train,y_test_all = train_test_split(x,y,test_size = 0.3,random_state=9)
    x_dev,x_test,y_dev,y_test = train_test_split(x_test_all,y_test_all,test_size=0.3,random_state=9)

    build_model(x_train,y_train,x_dev,y_dev)

它是如何工作的……

让我们从主函数开始。我们将调用get_data来获取我们的x预测属性和y响应属性。在get_data中，我们将利用make_classification数据集来生成用于随机森林方法的训练数据：

def get_data():
    """
    Make a sample classification dataset
    Returns : Independent variable y, dependent variable x
    """
    no_features = 30
    redundant_features = int(0.1*no_features)
    informative_features = int(0.6*no_features)
    repeated_features = int(0.1*no_features)
    x,y = make_classification(n_samples=500,n_features=no_features,flip_y=0.03,\
            n_informative = informative_features, n_redundant = redundant_features \
            ,n_repeated = repeated_features,random_state=7)
    return x,y

让我们来看一下传递给make_classification方法的参数。第一个参数是所需的实例数量。在这种情况下，我们需要 500 个实例。第二个参数是每个实例所需的属性数量。我们设置为 30。第三个参数flip_y，随机交换 3%的实例。这是为了在数据中引入噪声。接下来的参数是关于这 30 个特征中有多少个应该足够信息量，可以用于分类。我们指定 60%的特征，也就是 30 个中的 18 个，应该是有信息量的。接下来的参数是冗余特征。这些特征是通过信息特征的线性组合生成的，用来在特征之间引入相关性。最后，重复特征是从信息特征和冗余特征中随机抽取的重复特征。

让我们使用train_test_split将数据分为训练集和测试集。我们将预留 30%的数据用于测试：

    # Divide the data into Train, dev and test    
    x_train,x_test_all,y_train,y_test_all = train_test_split(x,y,test_size = 0.3,random_state=9)

再一次，我们将利用train_test_split将测试数据分成开发集和测试集：

    x_dev,x_test,y_dev,y_test = train_test_split(x_test_all,y_test_all,test_size=0.3,random_state=9)

我们将数据分为用于构建、评估和测试模型的三部分，然后开始构建我们的模型：

build_model(x_train,y_train,x_dev,y_dev)

在build_model中，我们将利用 scikit-learn 的SGDClassifier类来构建我们的随机梯度下降方法：

    estimator = SGDClassifier(n_iter=50,shuffle=True,loss="log", \
                learning_rate = "constant",eta0=0.0001,fit_intercept=True, penalty="none")

让我们来看一下我们使用的参数。第一个参数是我们希望遍历数据集的次数，以更新权重。在这里，我们设定为 50 次迭代。与感知机一样，在遍历完所有记录一次后，我们需要对输入记录进行洗牌，开始下一轮迭代。shuffle参数用于控制这个操作。shuffle的默认值是 true，这里我们加上它是为了说明。我们的损失函数是对数损失：我们希望进行逻辑回归，并且通过loss参数来指定这一点。我们的学习率 eta 是一个常数，我们将通过learning_rate参数来指定。我们将通过eta0参数来提供学习率的值。接着，我们将设置需要拟合截距，因为我们没有通过均值对数据进行中心化。最后，penalty参数控制所需的收缩类型。在我们的案例中，我们将设置不需要任何收缩，使用none字符串。

我们将通过调用fit函数来构建我们的模型，并用训练集和开发集来评估我们的模型：

 estimator.fit(x,y)
    train_predcited = estimator.predict(x)
    train_score = accuracy_score(y,train_predcited)
    dev_predicted = estimator.predict(x_dev)
    dev_score = accuracy_score(y_dev,dev_predicted)

    print 
    print "Training Accuracy = %0.2f Dev Accuracy = %0.2f"%(train_score,dev_score)

让我们看看我们的准确度得分：

还有更多内容……

对于 SGD 分类，可以应用正则化（L1、L2 或弹性网）。这个过程与回归相同，因此我们在这里不再重复。请参考前面的步骤。

在我们的示例中，学习率（eta）是常数，但这不一定是必要的。在每次迭代时，eta 值可以被减小。学习率参数learning_rate可以设置为一个最优的字符串或 invscaling。请参考以下 scikit 文档：

scikit-learn.org/stable/modules/sgd.html

该参数的设置如下：

estimator = SGDClassifier(n_iter=50,shuffle=True,loss="log", \
learning_rate = "invscaling",eta0=0.001,fit_intercept=True, penalty="none")

我们使用fit方法来构建我们的模型。如前所述，在大规模机器学习中，我们知道所有数据不会一次性提供给我们。当我们按批次接收数据时，需要使用partial_fit方法，而不是fit。使用fit方法会重新初始化权重，并且我们将失去上一批数据的所有训练信息。有关partial_fit的更多信息，请参考以下链接：

scikit-learn.org/stable/modules/generated/sklearn.linear_model.SGDClassifier.html#sklearn.linear_model.SGDClassifier.partial_fit

另请参见

• 在第七章中，使用岭回归进行收缩的食谱，机器学习 II

• 在第九章中，使用随机梯度下降进行回归的食谱，机器学习 III