Python机器学习实战：维数约简的6种常用方法

机器学习实战：这里没有艰深晦涩的数学理论，我们将用简单的案例和大量的示例代码，向大家介绍机器学习的核心概念。我们的目标是教会大家用Python构建机器学习模型，解决现实世界的难题。

为什么要约简维数？

当数据集包含大量特征，有的特征包含大量预测信息，有的仅包含少量信息或纯粹是噪音，很多特征之间也可能高度相关。维数约简的目的在于剔除噪音，只保留有意义的特征，这不仅使数据集更容易管理和理解，预测模型的准确性也会相应提高。

常用方法有哪些？

1. 相关系数矩阵
2. 主成分分析
3. 随机PCA
4. 因子分析
5. 线性判别分析
6. 核PCA

1. 相关系数矩阵

相关系数：度量数值变量同向/反向运动的关联程度，取值范围$[-1,1]$，1表示完全正相关，-1表示完全负相关，0表示线性无关。

Python实现：pd.corr(), np.corrcoef()用于计算相关系数矩阵，sns.heatmap()创建热力图可视化。

解读：如果两个特征高度相关，可以保留其中一个，尽可能保留相互独立的变量。

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import seaborn as sns

%matplotlib inline
plt.style.use("ggplot")

使用iris数据集。

from sklearn.datasets import load_iris

iris = load_iris()
iris_df = pd.DataFrame(iris.data, columns=iris.feature_names)
iris_df.head()

	sepal length (cm)	sepal width (cm)	petal length (cm)	petal width (cm)
0	5.1	3.5	1.4	0.2
1	4.9	3.0	1.4	0.2
2	4.7	3.2	1.3	0.2
3	4.6	3.1	1.5	0.2
4	5.0	3.6	1.4	0.2

计算相关系数矩阵。

cor_matrix = iris_df.corr()
np.round(cor_matrix, 2)

	sepal length (cm)	sepal width (cm)	petal length (cm)	petal width (cm)
sepal length (cm)	1.00	-0.12	0.87	0.82
sepal width (cm)	-0.12	1.00	-0.43	-0.37
petal length (cm)	0.87	-0.43	1.00	0.96
petal width (cm)	0.82	-0.37	0.96	1.00

创建热力图，可视化相关系数矩阵。

fig, ax = plt.subplots(figsize=(10, 8))
ax = sns.heatmap(
    cor_matrix,
    annot=True,  # 将数值添加到单元格
    vmin=-1,  # 根据数值范围调整颜色区间
    vmax=1, 
    center=0,
    linewidths=0.5,  # 控制单元格之间的间隙
    #mask=np.abs(cor_matrix) < 0.6,  # 仅显示相关系数大于0.6的单元格
    ax=ax
)

热力图显示，sepal length和petal length，petal width高度正相关，petal length和petal width高度相关。所以，可以考虑剔除petal length和petal width两个特征。

这个简单的例子仅用于说明如何用相关系数矩阵约简维数，考虑到原始数据只有4个特征，完全没有必要约简维数。

2. 主成分分析

主成分分析（PCA）是最常用的降维技术之一，它把大量相关变量转化为几组无关变量，这些无关变量称为“主成分”。

sklearn提供接口PCA类实现主成分分析。

拟合PCA模型前要先选择约简维数，即主成分个数，有两种常用方法：1. 参考相关文献的研究成果；2. 从相关系数矩阵获得启发。

from sklearn.decomposition import PCA

# 选择主成分个数，这里选择2
pca = PCA(n_components=2)
res = pca.fit_transform(iris_df)

print("original shape: ", iris_df.shape)
print("transformed shape: ", res.shape)

original shape:  (150, 4)
transformed shape:  (150, 2)

调用pca.fit_transform得到一个n维数组，每一列代表一个“主成分”，它们是原始特征的线性组合。

第一个主成分对原始数据的变异（方差）解释最大，第二个主成分与第一个主成分正交（无相关关系），它包含了剩余变异的大部分，第三个主成分与前两个主成分正交，并包含剩余变异的大部分，以此类推。

res[:5, :]

array([[-2.68412563,  0.31939725],
       [-2.71414169, -0.17700123],
       [-2.88899057, -0.14494943],
       [-2.74534286, -0.31829898],
       [-2.72871654,  0.32675451]])

用散点图查看主成分的关系，它们是线性无关的。

fig, ax = plt.subplots(figsize=(10, 7))
ax = sns.regplot(x=res[:, 0], y=res[:, 1], scatter=True, fit_reg=False, ax=ax)
ax.set_xlabel("First Main Component")
ax.set_ylabel("Second Main Component")

explained_variance_ratio_包含了所有主成分对原始数据方差的解释比率。降维不可避免地损失了原始数据的部分信息，降维的同时要尽可能多地保留信息。

var_ratio = pca.explained_variance_ratio_

for idx, val in enumerate(var_ratio, 1):
    print("principle component %d: %.2f%%" % (idx, val * 100))
print("total: %.2f%%" % np.sum(var_ratio * 100))

principle component 1: 92.46%
principle component 2: 5.31%
total: 97.77%

主成分个数等于2，就保留了原始数据集约98%的变异。

除了根据经验或相关系数矩阵外，可以用数据挖掘的方法选择主成分个数。拟合PCA模型时令主成分个数等于原始特征的数目，然后查看主成分的累积解释方差比率，选择令人满意的临界点对应的主成分个数。

# n_components=None, 主成分个数等于原始特征的数量
pca2 = PCA(n_components=None)
pca2.fit(iris_df)
evr = pca2.explained_variance_ratio_ * 100

for idx, val in enumerate(evr, 1):
    print("principle component %d: %.2f%%" % (idx, val))
print("total: %.2f%%" % np.sum(evr))

fig, ax = plt.subplots(figsize=(10, 7))
ax.plot(np.arange(1, len(evr) + 1), np.cumsum(evr), "r-o")
ax.set_title("Cumulative Explained Variance Ratio")
ax.set_xlabel("Number of components")
ax.set_ylabel("Explained Variance Ratio(%)")

principle component 1: 92.46%
principle component 2: 5.31%
principle component 3: 1.71%
principle component 4: 0.52%
total: 100.00%

也可以先选择最低的累计解释方差比率，让模型自动选择最优的主成分个数。

target = 0.95  # 最低累计解释方差比率
res = PCA(n_components=target).fit_transform(iris_df)
res[:5, :]

array([[-2.68412563,  0.31939725],
       [-2.71414169, -0.17700123],
       [-2.88899057, -0.14494943],
       [-2.74534286, -0.31829898],
       [-2.72871654,  0.32675451]])

3. 随机PCA

PCA类默认使用奇异值分解算法(SVD)，处理大数据时速度很慢，Randomized PCA采取基于随机SVD的算法，计算速度非常快，结果与经典SVD差异不大，特别适用于大型数据集。

先创建一个虚拟大型数据集，然后分别使用经典SVD和随机SVD算法，对比计算时间和最终结果。

from sklearn.datasets import make_classification
from sklearn.decomposition import PCA

X, y = make_classification(
    n_samples=10000,  # 10000个观测值
    n_features