【机器学习实战】-- Titanic 数据集（2）-- 感知机

1. 写在前面: 

本篇属于实战部分，更注重于算法在实际项目中的应用。如需对感知机算法本身有进一步的了解，可参考以下链接，在本人学习的过程中，起到了很大的帮助：

统计学习方法 李航

感知机原理小结 https://www.cnblogs.com/pinard/p/6042320.html

空间中任意一点到超平面距离的公式推导 https://www.cnblogs.com/yanganling/p/8007050.html

 

2. 数据集：

数据集地址：https://www.kaggle.com/c/titanic

Titanic数据集是Kaggle上参与人数最多的项目之一。数据本身简单小巧，适合初学者上手，深入了解比较各个机器学习算法。

数据集包含11个变量：PassengerID、Pclass、Name、Sex、Age、SibSp、Parch、Ticket、Fare、Cabin、Embarked，通过这些数据来预测乘客在Titanic事故中是否幸存下来。

 

3. 算法简介：

感知机属于分类模型，是一个古老而基础的模型，与支持向量机有一定程度的相似，同时也是神经网络的基础。

感知机属于线性模型，因此线性模型中常用的L1、L2正则化同样使用与感知机。

 

3.1 感知机模型：

由于不同材料中对多个多维数据的表达不尽相同，这里参考《统计学习方法》中李航老师的写法：

给定一个数据集：$T=\left \{ \left ( x_{1}, y_{1} \right ), \left ( x_{2}, y_{2} \right ), ..., \left ( x_{N}, y_{N} \right ) \right \}$，其中$x_{i}\in X\subseteq \bf{R^{n}}$，$y_{i} \in Y = \left \{+1, -1 \right \}$，$i = 1,2,...,N$。这代表数据集共有 N 对 实例，每个实例 $x_{i}$都是n维的。

 

从输入空间到输出空间的如下函数被称作感知机模型：

$f(x) = \rm{sign} \left( w \cdot x + b \right) $，其中sign是符号函数：$sign(x)= \begin{cases} +1& {x\geq0}\\ -1& {x< 0} \end{cases}$

 

3.2 感知机损失函数：

一般情况下，损失函数的选取是所有实例的预测值$f(x_{i})$与实际值$y_{i}$的差。而在感知机中，我们有一个特殊的损失函数，它更便于理解，也更便于后续的算法学习过程。我们将损失函数定义为：误分类点到超平面$S$的总距离。

首先，需要知道的是空间$\bf{R^{n}}$中任意一点$x_{0}$到超平面 $w \cdot x + b=0$的距离是：$ \frac{1}{\left \| w \right \|} \left | w \cdot x_{0} + b \right | $

其次，由于 $w \cdot x_{i} + b >0$时，$y_{i} = -1$，$w \cdot x_{i} + b <0$时，$y_{i} = 1$，因此，误分类点$x_{i}$到超平面$S$的距离是 $-\frac{1}{\left \| w \right \|} \left ( w \cdot x_{i} + b \right )$

由此，假设超平面$S$的误分类点集合为$M$，则所有误分类点到超平面$S$的总距离为：$-\frac{1}{\left \| w \right \|} \sum_{x_{i}\in M}y_{i} \left ( w \cdot x_{i} + b \right )$，感知机中只考虑函数间隔（即不考虑分母$\frac{1}{\left \| w \right \|}$）

最终，可以得到感知机的损失函数为：

$$ L\left ( w, b \right ) = -\sum_{x_{i} \in M} y_{i} \left ( w \cdot x_{i} + b \right ) $$

 3.3 感知机学习过程

感知机的学习过程即求解感知机损失函数的最小化问题：

$$ \min_{w,b} L(w,b) = -\sum_{x_{i} \in M} y_{i} (w \cdot x_{i} + b) $$

求解过程采用的是随机梯度下降 -- SGD(stochastic gradient descent)，损失函数的梯度为：

$$ \nabla_{w} L(w,b) = - \sum_{x_{i} \in M}y_{i}x_{i} $$

$$ \nabla_{b} L(w,b) = -\sum_{x_{i} \in M}y_{i} $$

每次选取随机的一个误分类点$(x_{i}, y_{i})$，对 $w$, $b$进行更新：

$$ w = w + \alpha y_{i}x_{i} $$

$$ b = b + \alpha y_{i} $$

其中 $\alpha$ 是学习率（learning rate)。

如果原原数据集是线性可分的，在合适的学习率下，经过迭代，损失函数$ L(w,b) $将减小至0。如果原数据是线性不可分的话，则需要设定其他的迭代终止条件。

3.4 感知机学习的对偶形式

感知机的对偶形式和原始形式相对应，但能提高求解速度。下面具体解释一下它是怎么提高求解速度的。

首先，假设$w$和$b$的初始取值都为0。对于误分类点$ (x_{i}, y_{i}) $，在学习过程中更新了$n_{i}$次，则 $w$, $b$ 关于$(x_{i}, y_{i})$的增量分别为 $n_{i}\alpha y_{i}x_{i}$ 和 $n_{i}\alpha y_{i}$。

接下来，考虑所有误分类点，并令 $\eta_{i} = n_{i}\alpha$，则最后学习到的 $w$ 和 $b$ 分别为：

$$ w = \sum_{i =1}^{N} \eta_{i}y_{i}x_{i} $$

$$ b = \sum_{i=1}^{N} \eta_{i}y_{i} $$

最后，有了以上结果，将$w$带回感知机模型 $f(x) =  \rm{sign}\left ( \sum_{j=1}^{N} \alpha_{j}y_{j}x_{j} \cdot x + b  \right ) $：

　　*) 每次迭代过程中，判断 $(x_{i}, y_{i})$ 是否为误分类点时，判断方程为：$y_{i}\left ( w \cdot x_{i} + b \right ) = y_{i}\left ( \sum_{j=1}^{N} \eta_{j}y_{j}x_{j} \cdot x_{i} + b \right )$

　　*) 如果判断方程 $y_{i}\left ( \sum_{j=1}^{N} \eta_{j}y_{j}x_{j} \cdot x_{i} + b \right ) \leq 0$，则更新 $\eta_{i} = \eta_{i} + \alpha$, $b = b + \alpha y_{i}$

对偶形式中能够提升求解速度的原因啊就在于，在对偶形式中，训练实例仅以内积的形式出现。预先将训练集中的实例间的内积计算出来，后续每次迭代过程中直接使用。而这个矩阵就是Gram矩阵：

$$ G=\left [ x_{i} \cdot x_{j} \right ]_{N\times N} $$

4. 实战： 

在本篇实战中，将在上一篇 【机器学习实战】-- Titanic 数据集（1）-- 朴素贝叶斯 的基础上，进一步介绍sklearn库在实战中的应用：

1. 使用Pipeline将数据清洗和模型预测结合起来，使代码更简洁（但同时也会带来麻烦）
2. 使用GridSearch选取最优参数
3. Perceptron的正则化

 首先，和之前一样，将数据进行导入：

import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
from sklearn.preprocessing import OneHotEncoder, OrdinalEncoder, StandardScaler
from sklearn.impute import SimpleImputer
from sklearn.pipeline import Pipeline, FeatureUnion
from sklearn.model_selection import cross_val_score, GridSearchCV, ParameterGrid, StratifiedKFold
from sklearn.linear_model import Perceptron
from sklearn.metrics import accuracy_score, precision_score, recall_score
from sklearn.base import TransformerMixin, BaseEstimator


class DataFrameSelector(BaseEstimator, TransformerMixin):
    def __init__(self, attribute_name):
        self.attribute_name = attribute_name

    def fit(self, x, y=None):
        return self

    def transform(self, x):
        return x[self.attribute_name].values


# Load data
data_train = pd.read_csv('train.csv')

train_x = data_train.drop('Survived', axis=1)
train_y = data_train['Survived']

 

 

4.1 数据清洗和模型预测结合起来

 

# Data cleaning
cat_attribs = ['Pclass', 'Sex', 'Embarked']
dis_attribs = ['SibSp', 'Parch']
con_attribs = ['Age', 'Fare']

# encoder: OneHotEncoder()、OrdinalEncoder()
cat_pipeline = Pipeline([
    ('selector', DataFrameSelector(cat_attribs)),
    ('imputer', SimpleImputer(strategy='most_frequent')),
    ('encoder', OneHotEncoder()),
])

dis_pipeline = Pipeline([
    ('selector', DataFrameSelector(dis_attribs)),
    ('scaler', StandardScaler()),
    ('imputer', SimpleImputer(strategy='most_frequent')),
])

con_pipeline = Pipeline([
    ('selector', DataFrameSelector(con_attribs)),
    ('scaler', StandardScaler()),
    ('imputer', SimpleImputer(strategy='mean')),
])

full_pipeline = FeatureUnion(
    transformer_list=[
        ('con_pipeline', con_pipeline),
        ('dis_pipeline', dis_pipeline),
        ('cat_pipeline', cat_pipeline),
    ]
)

# Applying the model without knowing train_x_cleaned
full_pipeline_estimator = Pipeline([
    ('clean', full_pipeline),
    ('predict', Perceptron()),
])

cv = StratifiedKFold(n_splits=5, shuffle=True, random_state=2)

param_grid = [
    {'predict__penalty': ['l2', 'l1', 'elasticnet', None],
     'predict__alpha': [0.002, 0.005, 0.0001],
     'predict__class_weight': ['balanced', None]},
]

grid_search = GridSearchCV(full_pipeline_estimator, param_grid, cv=cv, scoring='accuracy')
grid_search.fit(train_x, train_y)
predicted_y = grid_search.predict(train_x)

 

将数据清洗和模型预测结合起来，看似十分简洁，但其中隐藏着几个问题：

Q1. full_pipeline_estimator 无法直观地获取清洗后的数据，因为清洗只是整个pipeline中的一部分，直接使用 full_pipeline_estimator.fit_transform(train_x, train_y)甚至会报错，因为整个pipeline中最后一步Perception()并没有fit_transform()的用法。

 

Q2. 强行将 数据清洗 和 模型预测 结合在一起，实际上并没有对整个算法运行、参数调优起到什么作用。可以看到，在 GridSearchCV 的 param_grid 参数中，都是属于 ‘predict’ 这一步的参数。

 

Q3. 由于在数据清洗中，有涉及到缺失值的处理，无论是 ‘mean’ 还是 ‘most frequent’ 都是对全体数据而言的。而如果在 GridSearchCV 中 使用，程序将对每个fold中的数据单独进行数据清理，将导致每个fold中对缺失值的处理不一致。

 

对于第1，第2个问题，我们可以参考sklearn的官方文档获得相应解答：

 

“The purpose of the pipeline is to assemble several steps that can be cross-validated together while setting different parameters. For this, it enables setting parameters of the various steps using their names and the parameter name separated by a ‘__’. ”

Pipeline的目的，是将若干个可以交叉验证的步骤组合在一起，并且同时可以设置不同的参数。为此，它能够通过 “步骤名_参数名” 的形式来设置参数。

 

A1：因此，对于问题1，可以通过以来代码来获取清洗后的数据：

 train_x_cleaned = grid_search.best_estimator_.named_steps['clean'].fit_transform(train_x) 

 

A2：对于问题2，从官方文档中也获得了解答：无需设置参数的步骤，不需要放入pipeline中。那我们可以将什么放入piepline中呢，比如数据清洗后的降维，如PCA等。

 

A3：第三个问题则是比较严重的问题。

首先，我们查看整个 grid_search 的最优解：

 

In[2]:
grid_search.best_params_
Out[2]: 
{'predict__alpha': 0.0001,
 'predict__class_weight': None,
 'predict__penalty': 'l2'}

 

In[21]:
print('Result of full pipeline estimator:')
print(accuracy_score(train_y, predicted_y))
print(precision_score(train_y, predicted_y))
print(recall_score(train_y, predicted_y))
Out[21]:
Result of full pipeline estimator:
0.7833894500561167
0.752542372881356
0.6491228070175439

 

接下来，对每个 fold中的最优解进行查看，查看每个fold中的最优解和整个grid_search的最优解是否一致。在这一过程中，我们希望找到的最优解是针对这个fold中的所有数据的，所以不再需要 cross-validation。有两种方法可以达成这一效果：

（1）巧妙地使用GridSearchCV中的关键字cv，使用 cv2= [(slice(None), slice(None))]

In[48]:
for train_index, test_index in cv.split(train_x, train_y):
    X_train, X_test = train_x.loc[train_index], train_x.loc[test_index]
    y_train, y_test = train_y.loc[train_index], train_y.loc[test_index]
    cv2 = [(slice(None), slice(None))]
    grid_search = GridSearchCV(full_pipeline_estimator, param_grid, cv=cv2, scoring='accuracy')
    grid_search.fit(X_train, y_train)
    predicted_y = grid_search.predict(X_train)
    print(accuracy_score(y_train, predicted_y))
    print(grid_search.best_params_)
    
Out[48]:
0.7907303370786517
{'predict__alpha': 0.0001, 'predict__class_weight': None, 'predict__penalty': 'l1'}
0.7629733520336606
{'predict__alpha': 0.0001, 'predict__class_weight': None, 'predict__penalty': 'l1'}
0.7769985974754559
{'predict__alpha': 0.002, 'predict__class_weight': 'balanced', 'predict__penalty': None}
0.7980364656381487
{'predict__alpha': 0.0001, 'predict__class_weight': None, 'predict__penalty': 'l1'}
0.7699859747545582
{'predict__alpha': 0.0001, 'predict__class_weight': None, 'predict__penalty': 'l2'}

（2）使用 ParameterGrid，遍历其中的参数：

In[49]:
param_grid_ = [
    {'penalty': ['l2', 'l1', 'elasticnet', None],
     'alpha': [0.002, 0.005, 0.0001],
     'class_weight': ['balanced', None]},
]
for train_index, test_index in cv.split(train_x, train_y):
    X_train, X_test = train_x.loc[train_index], train_x.loc[test_index]
    y_train, y_test = train_y.loc[train_index], train_y.loc[test_index]
    best_score = 0
    best_parameter = []
    for g in ParameterGrid(param_grid_):
        full_pipeline_estimator.named_steps['predict'].set_params(**g)
        full_pipeline_estimator.fit(X_train, y_train)
        predicted_y = full_pipeline_estimator.predict(X_train)
        score = accuracy_score(y_train, predicted_y)
        if score > best_score:
            best_score = score
            best_parameter = g
    print(best_score)
    print(best_parameter)
    
Out[49]:
0.7907303370786517
{'alpha': 0.0001, 'class_weight': None, 'penalty': 'l1'}
0.7629733520336606
{'alpha': 0.0001, 'class_weight': None, 'penalty': 'l1'}
0.7769985974754559
{'alpha': 0.002, 'class_weight': 'balanced', 'penalty': None}
0.7980364656381487
{'alpha': 0.0001, 'class_weight': None, 'penalty': 'l1'}
0.7699859747545582
{'alpha': 0.0001, 'class_weight': None, 'penalty': 'l2'}

可以看到两者的效果是一致的。并且grid_search全局的最优解，只在最后一个fold中是最优解。这一方面可能是由于不同参数组合间的结果差距很小；另一方面，对数据预处理的不同，导致它们面对的并不是完全相同的训练数据

 

4.2 数据清洗和模型预测分开

在上一节中，我们看到了将数据清洗和模型预测合在一起的弊端，那么我们再将其分开，看看结果：

# Data cleaning
cat_attribs = ['Pclass', 'Sex', 'Embarked']
dis_attribs = ['SibSp', 'Parch']
con_attribs = ['Age', 'Fare']

# encoder: OneHotEncoder()、OrdinalEncoder()
cat_pipeline = Pipeline([
    ('selector', DataFrameSelector(cat_attribs)),
    ('imputer', SimpleImputer(strategy='most_frequent')),
    ('encoder', OneHotEncoder()),
])

dis_pipeline = Pipeline([
    ('selector', DataFrameSelector(dis_attribs)),
    ('scaler', StandardScaler()),
    ('imputer', SimpleImputer(strategy='most_frequent')),
])

con_pipeline = Pipeline([
    ('selector', DataFrameSelector(con_attribs)),
    ('scaler', StandardScaler()),
    ('imputer', SimpleImputer(strategy='mean')),
])

full_pipeline = FeatureUnion(
    transformer_list=[
        ('con_pipeline', con_pipeline),
        ('dis_pipeline', dis_pipeline),
        ('cat_pipeline', cat_pipeline),
    ]
)

# Applying the model with train_x_cleaned
train_x_cleaned = full_pipeline.fit_transform(train_x)

clf = Perceptron()
# cv = StratifiedKFold(n_splits=5, shuffle=True, random_state=2)
# cross_val_score(clf, train_x_cleaned, train_y, cv=cv)

param_grid2 = [
    {'penalty': ['l2', 'l1', 'elasticnet', None],
     'alpha': [0.002, 0.005, 0.0001],
     'class_weight': ['balanced', None]},
]

grid_search2 = GridSearchCV(clf, param_grid2, cv=cv, scoring='accuracy')
grid_search2.fit(train_x_cleaned, train_y)

predicted_y2 = grid_search2.predict(train_x_cleaned)

查看整个 grid_search 的最优解：

In[8]: grid_search2.best_params_
Out[8]: {'alpha': 0.0001, 'class_weight': None, 'penalty': 'l1'}

In[22]:
print('Result of full pipeline + separate estimator:')
print(accuracy_score(train_y, predicted_y2))
print(precision_score(train_y, predicted_y2))
print(recall_score(train_y, predicted_y2))
Out[22]:
Result of full pipeline + separate estimator:
0.7710437710437711
0.6938202247191011
0.7222222222222222

 

同样，对每个 Fold 中的最优解进行查看：

In[55]:
for train_index, test_index in cv.split(train_x_cleaned, train_y):
    X_train, X_test = train_x_cleaned[train_index], train_x_cleaned[test_index]
    y_train, y_test = train_y[train_index], train_y[test_index]
    cv2 = [(slice(None), slice(None))]
    # cv = ShuffleSplit(test_size=0.2, train_size=1.0, n_splits=1, random_state=0)
    grid_search2 = GridSearchCV(clf, param_grid2, cv=cv2, scoring='accuracy')
    grid_search2.fit(X_train, y_train)
    predicted_y = grid_search2.predict(X_train)
    print(accuracy_score(y_train, predicted_y))
    print(grid_search2.best_params_)
    
Out[55]:
0.8061797752808989
{'alpha': 0.0001, 'class_weight': None, 'penalty': 'l1'}
0.7924263674614306
{'alpha': 0.005, 'class_weight': None, 'penalty': 'l1'}
0.7980364656381487
{'alpha': 0.002, 'class_weight': 'balanced', 'penalty': 'l2'}
0.791023842917251
{'alpha': 0.002, 'class_weight': None, 'penalty': 'l2'}
0.7741935483870968
{'alpha': 0.0001, 'class_weight': None, 'penalty': 'l1'}

可以看到，整个 grid_search 的最优解，在2个fold的数据中是最优解，占比较“数据清洗+模型预测”的方式有所提升，这一定程度上是由于数据预处理的相同。

 

4.3 Perceptron的正则化

在上两节中，可以看到 Perception的参数中有 'l1'， 'l2'。它们都是正则化的一种。同时可以通过 coef_ 参数发现，‘l1’ 正则化，会有许多系数是0，而 ‘l2’ 正则化则没有。这体现了 ‘l1’ 正则化的一个特点：可以在约束参数的同时，对数据进行降维。

 

4.4 结果分析

我们将4.2节中的方法得到的最优参数 " alpha = 0.0001, class_weight = None, penalty = 'l1' "，用于预测集，并将结果上传kaggle，结果如下：

	训练集 accuracy	训练集 precision	训练集 recall	预测集 accuracy（需上传kaggle获取结果）
朴素贝叶斯最优解	0.790	0.731	0.716	0.756
感知机	0.771	0.694	0.722	0.722