针对科学数据处理的统计学习教程（scikit-learn教程2）

http://www.cnblogs.com/taceywong/p/4570155.html

统计学习

：

随着科学实验数据的迅速增长，机器学习成了一种越来越重要的技术。问题从构建一个预测函数将不同的观察数据联系起来，到将观测数据分类，或者从未标记数据中学习到一些结构。

本教程将探索机器学习中统计推理的统计学习的使用：将手中的数据做出结论

Scikit-learn

 是一个紧密结合Python科学计算库(Numpy、Scipy、matplotlib)，集成经典机器学习算法的Python模块。

一、统计学习：scikit-learn中的设置与评估函数对象

（1）数据集

scikit-learn 从二维数组描述的数据中学习信息。他们可以被理解成多维观测数据的列表。如（n,m），n表示样例轴，y表示特征轴。

使用scikit-learn装载一个简单的样例：iris数据集

>>from sklearn import datasets
>>iris = datasets.load_iris()
>>data = iris.data
>>data.shape
(150, 4)

它有150个iris观测数据构成，每一个样例被四个特征所描述：他们的萼片、花瓣长度、花瓣宽度，具体的信息可以通过iris》DESCR查看。

当数据初始时不是(

n样例，n特征

)样式时，需要将其预处理以被scikit-learn使用。

通过数字数据集讲述数据变形

数字数据集由1797个8x8手写数字图片组成

>>>digits = datasets.load_digits()
>>>digits.images.shape
(1797, 8, 8)
>>> import pylab as pl
>>>pl.imshow(digits.images[-1], cmap=pl.cm.gray_r)
<matplotlib.image.AxesImage object at ...>

在scikit-learn中使用这个数据集，我们需要将其每一个8x8图片转换成长64的特征向量

python
>>>data = digits.images.reshape((digits.images.shape[0],-1))

(2)估计函数对象

拟合数据

：scikit-learn实现的主要API是估计函数。估计函数是用以从数据中学习的对象。它可能是分类、回归、聚类算法，或者提取过滤数据特征的转换器。

一个估计函数带有一个

fit

方法，以dataset作为参数（一般是个二维数组）

>>>estimator.fit(data)

估计函数对象的参数

：每一个估测器对象在实例化或者修改其相应的属性，其参数都会被设置。

>>>estimator = Estimator(param1=1, param2=2)
>>>estimator.param1
1

估测后的参数

：

>>>estimator.estimated_param_

二、有监督学习：从高维观察数据预测输出变量

有监督学习解决的问题

有监督学习主要是学习将两个数据集联系起来：观察数据x和我们要尝试预测的外置变量y，y通常也被称作目标、标签。多数情况下，y是一个和n个观测样例对应的一维数组。

scikit-learn中实现的所有有监督学习评估对象，都有fit(X,Y)方法来拟合模型，predict(X)方法根据未加标签的观测数据X

返回预测的标签y。

词汇：分类和回归

如果预测任务是将观测数据分类到一个有限的类别集中，换句话说，给观测对象命名，那么这个任务被称作分类任务。另一方面，如果任务的目标是预测测目标是一个连续性变量，那么这个任务成为回归任务。

用scikit-learn解决分类问题时，y是一个整数或字符串组成的向量

注意：查看[]快速了解用scikit-learn解决机器学习问题过程中的基础词汇。

（1）近邻和高维灾难

iris分类

：

iris分类是根据花瓣、萼片长度、萼片宽度来识别三种不同类型的iris的分类任务:

>> import numpy as np
>> from sklearn import datasets
>> iris = datasets.load_iris()
>> iris_X = iris.data
>> iris_y = iris.target
>> np.unique(iris_y)
array([0, 1, 2])

最近邻分类器

：

近邻也许是最简的分类器：得到一个新的观测数据X-test，从训练集的观测数据中寻找特征最相近的向量。（【】）

训练集和测试集

：

当尝试任何学习算法的时候，评估一个学习算法 的预测精度是很重要的。所以在做机器学习相关的问题的时候，通常将数据集分成训练集和测试集。

KNN(最近邻)分类示例：

# Split iris data in train and test data
# A random permutation, to split the data randomly
np.random.seed(0)
indices = np.random.permutation(len(iris_X))
iris_X_train = iris_X[indices[:-10]]
iris_y_train = iris_y[indices[:-10]]
iris_X_test  = iris_X[indices[-10:]]
iris_y_test  = iris_y[indices[-10:]]
# Create and fit a nearest-neighbor classifier
from sklearn.neighbors import KNeighborsClassifier
knn = KNeighborsClassifier()
knn.fit(iris_X_train, iris_y_train)
knn.predict(iris_X_test)
iris_y_test

高维灾难：

对于一个有效的学习算法，你需要最近n个点之间的距离d（依赖于具体的问题）。在一维空间中，需要平局n1/d各点，在上文中提到的K-NN例子中，如果数据只是有一个0-1之间的特征和n个训练观测数据所表述的画，那么新数据将不会超过1/n。因此，最近邻决策规则非常高效，因为与类间特征变化的范围相比，1/n小的多。

如果特征数是P，你就需要n 1/d^p个点。也就是说，如果我们在一维度情况下需要10个点，在P维度情况下需要10^p个点。当P变大的时候，为获得一个好的预测函数需要的点数将急剧增长。

这被称为高维灾难（指数级增长），也是机器学习领域的一个核心问题。

（2）线性模型：从回归到稀疏性

Diabets数据集（糖尿病数据集）
糖尿病数据集包含442个患者的10个生理特征（年龄，性别、体重、血压）和一年以后疾病级数指标。

diabetes = datasets.load_diabetes()

diabetes_X_train = diabetes.data[:-20]

diabetes_X_test = diabetes.data[-20:]

diabetes_y_train = diabetes.target[:-20]

diabetes_y_test = diabetes.target[-20:]

手上的任务是从生理特征预测疾病级数

线性回归：

【线性回归】的最简单形式给数据集拟合一个线性模型，主要是通过调整一系列的参以使得模型的残差平方和尽量小。

线性模型：y = βX+b
X:数据
y：目标变量
β：回归系数
b:观测噪声（bias，偏差）

from sklearn import linear_model
regr = linear_model.LinearRegression()
regr.fit(diabetes_X_train, diabetes_y_train)
print(regr.coef_)
# The mean square error
np.mean((regr.predict(diabetes_X_test)-diabetes_y_test)**2)
# Explained variance score: 1 is perfect prediction
# and 0 means that there is no linear relationship
# between X and Y.
regr.score(diabetes_X_test, diabetes_y_test)

收缩（Shrinkage）：

如果每一维的数据点很少，噪声将会造成很大的偏差影响：

X = np.c_[ .5, 1].T
y = [.5, 1]
test = np.c_[ 0, 2].T
regr = linear_model.LinearRegression()
import pylab as pl
pl.figure()
np.random.seed(0)
for _ in range(6):
this_X = .1*np.random.normal(size=(2, 1)) + X
regr.fit(this_X, y)
pl.plot(test, regr.predict(test))
pl.scatter(this_X, y, s=3)

高维统计学习的一个解决方案是将回归系数缩小到0：观测数据中随机选择的两个数据集近似不相关。这被称为岭回归（Ridge Regression）：

regr = linear_model.Ridge(alpha=.1)
pl.figure()
np.random.seed(0)
for _ in range(6):
this_X = .1*np.random.normal(size=(2, 1)) + X
regr.fit(this_X, y)
pl.plot(test, regr.predict(test))
pl.scatter(this_X, y, s=3)

这是一个偏差/方差（bias/variance）的权衡：岭α参数越大，偏差（bias)越大，方差（variance）越小

我们可以选择α以最小化排除错误，这里使用糖尿病数据集而不是人为制造的数据：

alphas = np.logspace(-4, -1, 6)
from __future__ import print_function
print([regr.set_params(alpha=alpha
).fit(diabetes_X_train, diabetes_y_train,
).score(diabetes_X_test, diabetes_y_test) for alpha in alphas])

【注意】扑捉拟合参数的噪声使得模型不能推广到新的数据被称为过拟合。岭回归造成的偏差被称为正则化（归整化，regularization）

稀疏性：

只拟合特征1和特征2：

【注意】整个糖尿病数据包含11维数据（

10个特征维，一个目标变量

），很难对这样的数据直观地表现出来，但是记住那是一个很空的空间也许是有用的。

我们可以看到，尽管特征2在整个模型中占据很大的系数，但是和特征1相比，对结果y造成的影响很小。

为了提升问题的状况（考虑到高维灾难），只选择信息含量较大的（对结果y造成的影响较大的）的特征，不选择信息含量较小的特征会很有趣，如把特征2的系数调到0.岭回归将会减少信息含量较小的系数的值，而不是把它们设置为0.另一种抑制措施——Lasso（最小绝对收缩和选择算子）可以使得一些参数为0.这些方法被称作稀疏方法。系数操作可以看作是奥卡姆的剃刀：模型越简单越好。

regr = linear_model.Lasso()
scores = [regr.set_params(alpha=alpha
).fit(diabetes_X_train, diabetes_y_train
).score(diabetes_X_test, diabetes_y_test)
for alpha in alphas]
best_alpha = alphas[scores.index(max(scores))]
regr.alpha = best_alpha
regr.fit(diabetes_X_train, diabetes_y_train)
print(regr.coef_)

针对相同问题的不同算法：

不同的算法可以被用来解决相同的数学问题。例如scikit-learn中的Lasso对象使用coordinate decent方法解决lasso回归问题，在大数据集上是很有效的。然而，scikit-learn也使用LARS算法提供了LassoLars对象，对于处理权重向量非常稀疏的数据非常有效（数据的观测实例非常少）。
分类：

对于分类问题，比如iris标定任务，线性回归不是正确的方法。因为它会给数据得出大量远离决策边界的权重。一个线性方法是你和一个sigmoid函数或者logistic函数：

logistic = linear_model.LogisticRegression(C=1e5)
logistic.fit(iris_X_train, iris_y_train)

这就是有名的logistic回归。
多分类：

如果你有多个类别需要预测，一个可行的方法是 “一对多”分类，接着根据投票决定最终的决策。

通过Logistic回归进行收缩和稀疏：

在LogisticRegression对象中C参数控制着正则化的数量：C越大，正则化数目越少。penalty= "12" 提供收缩（非稀疏化系数），penalty="11"提供稀疏化。

练习：

尝试使用近邻算法和线性模型对数字数据集进行分类。留出最后的10%作为测试集用来测试预测的精确度。

from sklearn import datasets, neighbors, linear_model
digits = datasets.load_digits()
X_digits = digits.data
y_digits = digits.target

【完整代码】

from sklearn import datasets, neighbors, linear_model
digits = datasets.load_digits()
X_digits = digits.data
y_digits = digits.target
n_samples = len(X_digits)
X_train = X_digits[:.9 * n_samples]
y_train = y_digits[:.9 * n_samples]
X_test = X_digits[.9 * n_samples:]
y_test = y_digits[.9 * n_samples:]
knn = neighbors.KNeighborsClassifier()
logistic = linear_model.LogisticRegression()
print('KNN score: %f' % knn.fit(X_train, y_train).score(X_test, y_test))
print('LogisticRegression score: %f'
% logistic.fit(X_train, y_train).score(X_test, y_test))

（3）支持向量机（SVMs）

线性SVNs：

支持向量机属于判别模型家族：它们尝试寻找样例的一个组合，构建一个两类之间的最大边缘平面。通过C参数进行正则化：一个较小的C意味着边缘是通过分割线周围的所有观测样例进行计算得到的（更规整化，正则化）；一个较大的C意味着边缘是通过邻近分割线的观测样例计算得到的（更少的规整化，正则化）：

非正则化SVN:

正则化 SVM(默认)：

样例：Plot different SVM分类器 iris数据集

SVMs能够被用于回归——SVR（支持向量回归）—用于分类——SVC(支持向量分类)

from sklearn import svm

svc = svm.SVC(kernel='linear')

svc.fit(iris_X_train, iris_y_train)
【警告】：规格化数据

对于大多数的估测模型，包括SVMs，处理好单位标准偏差对于获得一个好的预测是很重要的。

使用核函数：

在特征空间中类别不经常是线性可分的。解决方案是构建一个非线性但能用多项式代替的决策函数。这要通过核技巧实现：使用核可以被看作通过设置核在观测样例上创建决策力量。

线性核：

多项式核：

径向基函数核（RBF,Radial Basis Function）：

svc = svm.SVC(kernel='rbf')

交互式样例：

参照SVM GUI，下载svm_gui.py;通过鼠标左右键设置两类数据点，拟合模型并改变参数和数据。

练习：

尝试使用SVMs根据iris数据集前两个特征将其分成两类。留出每一类的10%作为测试样例。
【警告】数据集中的数据是按照分类顺序排列的，不要留出最后的10%作为测试样例，要不然你只能测试一种类别。（获取训练集和测试集是注意要进行混淆）
提示：你可以在一个网格上使用decision_function方法获得直观的呈现。

iris = datasets.load_iris()
X = iris.data
y = iris.target
X = X[y != 0, :2]
y = y[y != 0]

完整代码：

"""
================================
SVM Exercise
================================
A tutorial exercise for using different SVM kernels.
This exercise is used in the :ref:`using_kernels_tut` part of the
:ref:`supervised_learning_tut` section of the :ref:`stat_learn_tut_index`.
"""
print(__doc__)
import numpy as np
import matplotlib.pyplot as plt
from sklearn import datasets, svm
iris = datasets.load_iris()
X = iris.data
y = iris.target
X = X[y != 0, :2]
y = y[y != 0]
n_sample = len(X)
np.random.seed(0)
order = np.random.permutation(n_sample)
X = X[order]
y = y[order].astype(np.float)
X_train = X[:.9 * n_sample]
y_train = y[:.9 * n_sample]
X_test = X[.9 * n_sample:]
y_test = y[.9 * n_sample:]
# fit the model
for fig_num, kernel in enumerate(('linear', 'rbf', 'poly')):
clf = svm.SVC(kernel=kernel, gamma=10)
clf.fit(X_train, y_train)
plt.figure(fig_num)
plt.clf()
plt.scatter(X[:, 0], X[:, 1], c=y, zorder=10, cmap=plt.cm.Paired)
# Circle out the test data
plt.scatter(X_test[:, 0], X_test[:, 1], s=80, facecolors='none', zorder=10)
plt.axis('tight')
x_min = X[:, 0].min()
x_max = X[:, 0].max()
y_min = X[:, 1].min()
y_max = X[:, 1].max()
XX, YY = np.mgrid[x_min:x_max:200j, y_min:y_max:200j]
Z = clf.decision_function(np.c_[XX.ravel(), YY.ravel()])
# Put the result into a color plot
Z = Z.reshape(XX.shape)
plt.pcolormesh(XX, YY, Z > 0, cmap=plt.cm.Paired)
plt.contour(XX, YY, Z, colors=['k', 'k', 'k'], linestyles=['--', '-', '--'],
levels=[-.5, 0, .5])
plt.title(kernel)
plt.show()

三、模型选择：选择模型和他们的参数

（1）分数，和交叉验证分数

众所周知，每一个模型会得出一个score方法用于裁决模型在新的数据上拟合的质量。其值越大越好。

from sklearn import datasets, svm
digits = datasets.load_digits()
X_digits = digits.data
y_digits = digits.target
svc = svm.SVC(C=1, kernel='linear')
svc.fit(X_digits[:-100], y_digits[:-100]).score(X_digits[-100:], y_digits[-100:])

为了获得一个更好的预测精确度度量，我们可以把我们使用的数据折叠交错地分成训练集和测试集：

import numpy as np
X_folds = np.array_split(X_digits, 3)
y_folds = np.array_split(y_digits, 3)
scores = list()
for k in range(3):
# We use 'list' to copy, in order to 'pop' later on
X_train = list(X_folds)
X_test  = X_train.pop(k)
X_train = np.concatenate(X_train)
y_train = list(y_folds)
y_test  = y_train.pop(k)
y_train = np.concatenate(y_train)
scores.append(svc.fit(X_train, y_train).score(X_test, y_test))
print(scores)

这被称为KFold交叉验证

（2）交叉验证生成器

上面将数据划分为训练集和测试集的代码写起来很是沉闷乏味。scikit-learn为此自带了交叉验证生成器以生成目录列表：

from sklearn import cross_validation
k_fold = cross_validation.KFold(n=6, n_folds=3)
for train_indices, test_indices in k_fold:
print('Train: %s | test: %s' % (train_indices, test_indices))

接着交叉验证就可以很容易实现了：

kfold = cross_validation.KFold(len(X_digits), n_folds=3)
[svc.fit(X_digits[train], y_digits[train]).score(X_digits[test], y_digits[test])
for train, test in kfold]

为了计算一个模型的score，scikit-learn自带了一个帮助函数：

cross_validation.cross_val_score(svc, X_digits, y_digits, cv=kfold, n_jobs=-1)

n_jobs=-1

意味着将计算任务分派个计算机的所有CPU.

交叉验证生成器：

KFold(n,k)

 交叉分割，K-1上进行训练，生于数据样例用于测试

StratifiedKFold(y,K)

 保存每一个fold的类比率/标签分布

leaveOneOut(n)

 至预留一个观测样例

leaveOneLabelOut(labels)

 采用一个标签数组把观测样例分组

练习：

使用digits数据集，绘制使用线性核的SVC进行交叉验证的分数（使用对数坐标轴，1——10）

import numpy as np
from sklearn import cross_validation, datasets, svm
digits = datasets.load_digits()
X = digits.data
y = digits.target
svc = svm.SVC(kernel='linear')
C_s = np.logspace(-10, 0, 10)

完整代码：

（3）网格搜索和交叉验证模型

网格搜索：

scikit-learn提供一个对象，他得到数据可以在采用一个参数的模型拟合过程中选择使得交叉验证分数最高的参数。该对象的构造函数需要一个模型作为参数：

from sklearn.grid_search import GridSearchCV
Cs = np.logspace(-6, -1, 10)
clf = GridSearchCV(estimator=svc, param_grid=dict(C=Cs),
n_jobs=-1)
clf.fit(X_digits[:1000], y_digits[:1000])
clf.best_score_
clf.best_estimator_.C
# Prediction performance on test set is not as good as on train set
clf.score(X_digits[1000:], y_digits[1000:])

默认情况下，

GridSearchCV

使用

3-fold

交叉验证。然而，当他探测到是一个分类器而不是回归量，将会采用分层的

3-fold

。
嵌套 交叉验证

cross_validation.cross_val_score(clf, X_digits, y_digits)

两个交叉验证循环是并行执行的：一个

GridSearchCV

模型设置

gamma

,另一个使用

cross_val_score

 度量模型的预测表现。结果分数是在新数据预测分数的无偏差估测。

【警告】你不能在并行计算时嵌套对象（

n_jobs

不同于1）

交叉验证估测：

在算法by算法的基础上使用交叉验证去设置参数更高效。这也是为什么对于一个特定的模型/估测器引入

Cross-validation

:评估估测器表现模型去自动的通过交叉验证设置参数。

from sklearn import linear_model, datasets
lasso = linear_model.LassoCV()
diabetes = datasets.load_diabetes()
X_diabetes = diabetes.data
y_diabetes = diabetes.target
lasso.fit(X_diabetes, y_diabetes)
# The estimator chose automatically its lambda:
lasso.alpha_

这些模型的称呼和他们的对应模型很相似，只是在他们模型名字的后面加上了'

CV

'.

练习：

使用糖尿病数据集，寻找最佳的正则化参数

α

附加：你对选择的α值信任度有多高？

from sklearn import cross_validation, datasets, linear_model
diabetes = datasets.load_diabetes()
X = diabetes.data[:150]
y = diabetes.target[:150]
lasso = linear_model.Lasso()
alphas = np.logspace(-4, -.5, 30)

完整代码：

```python

```

四、无监督学习：寻找数据的代表

（1）聚类：将观测样例聚集到一起

聚类解决的问题：

比如对于iris数据集，如果我们知道我们知道有三种iris，但是我们没有标签标定他们：我们可以尝试聚类任务：将观测样例分成分离的族群中，这些族群可以被称为簇。
K-mean聚类（K均值聚类）

注意存在很多不同的聚类标准和关联算法。最简的聚类算法是——K均值（K-means）

from sklearn import cluster, datasets
iris = datasets.load_iris()
X_iris = iris.data
y_iris = iris.target
k_means = cluster.KMeans(n_clusters=3)
k_means.fit(X_iris)
print(k_means.labels_[::10])
print(y_iris[::10])

注意：没有绝对的保证能够恢复真实的分类。首先，尽管scikit-learn使用很多技巧来缓和问题的难度，但选择簇的个数还是是很困难的，初始状态下算法是很敏感的，可能会陷入局部最小。

不好的初始状态：

8个簇：

真实情况：

不要“过解释”聚类结果

应用实例：矢量化

K-means和一般的聚类，可以看作是选择少量的示例压缩信息的方式。这个问题被称之为矢量化。例如，这可以被用于分离一个图像：

import scipy as sp
try:
lena = sp.lena()
except AttributeError:
from scipy import misc
lena = misc.lena()
X = lena.reshape((-1, 1)) # We need an (n_sample, n_feature) array
k_means = cluster.KMeans(n_clusters=5, n_init=1)
k_means.fit(X)
values = k_means.cluster_centers_.squeeze()
labels = k_means.labels_
lena_compressed = np.choose(labels, values)
lena_compressed.shape = lena.shape

原始图像：

K-means矢量化：

等段：（Equal bins）

图像直方图：
分层凝聚聚类：Ward

分层聚类方法是一种针对构建一个簇的分层的簇分析。通常它的实现方式有以下两种：
凝聚：自下而上的方法：每一个观测样例开始于他自己的簇，以一种最小连接标准迭代合并。这种方法在观测样例较少的情况下非常有效（有趣）。当簇的数量变大时，计算效率比K-means高的多。
分裂：自上而下的方法：所有的观测样例开始于同一个簇。迭代的进行分层。对于预计簇很多的情况，这种方法既慢（由于所有的观测样例作为一个簇开始的，是递归进行分离的）又有统计学行的病态。

连同-驱使聚类（Conectivity-constrained clustering）

使用凝聚聚类，通过一个连通图可以指定某些样例能被聚集在一起。scikit-learn中的图通过邻接矩阵来表示，且通常是一个稀疏矩阵。例如，在聚类一张图片时检索连通区域（有时也被称作连同单元、部件）：

from sklearn.feature_extraction.image import grid_to_graph
from sklearn.cluster import AgglomerativeClustering
###############################################################################
# Generate data
lena = sp.misc.lena()
# Downsample the image by a factor of 4
lena = lena[::2, ::2] + lena[1::2, ::2] + lena[::2, 1::2] + lena[1::2, 1::2]
X = np.reshape(lena, (-1, 1))
###############################################################################
# Define the structure A of the data. Pixels connected to their neighbors.
connectivity = grid_to_graph(*lena.shape)
###############################################################################
# Compute clustering
print("Compute structured hierarchical clustering...")
st = time.time()
n_clusters = 15  # number of regions
ward = AgglomerativeClustering(n_clusters=n_clusters,
linkage='ward', connectivity=connectivity).fit(X)
label = np.reshape(ward.labels_, lena.shape)
print("Elapsed time: ", time.time() - st)
print("Number of pixels: ", label.size)
print("Number of clusters: ", np.unique(label).size)

特征凝聚：

我们已经知道稀疏性可以缓和高维灾难。i.e相对于特征数量观测样例数量不足的情况。另一种方法是合并相似的特征：特征凝聚。这种方法通过在特征方向上进行聚类实现。在特征方向上聚类也可以理解为聚合转置的数据。

digits = datasets.load_digits()
images = digits.images
X = np.reshape(images, (len(images), -1))
connectivity = grid_to_graph(*images[0].shape)
agglo = cluster.FeatureAgglomeration(connectivity=connectivity,
n_clusters=32)
agglo.fit(X)
X_reduced = agglo.transform(X)
X_approx = agglo.inverse_transform(X_reduced)
images_approx = np.reshape(X_approx, images.shape)

transeform

 和

invers_transeform

方法

有些模型带有转置方法。例如用来降低数据集的维度

（2）分解：从一个信号到成分和加载

成分及其加载：

如果X是我们的多变量数据，那么我们要要尝试解决的问题就是在不同的观测样例上复写写它：我们想要学习加载L和其它一系列的成分C，如X = LC。存在不同的标准和条件去选择成分。
主成分分析：PCA

主成分分析（PCA）选择在信号上解释极大方差的连续成分。

上面观测样例的点分布在一个方向上是非常平坦的：三个特征单变量的一个甚至可以有其他两个准确的计算出来。PCA用来发现数据在哪个方向上是不平坦的。

当被用来转换数据的时候，PCA可以通过投射到一个主子空间来降低数据的维度。：

# Create a signal with only 2 useful dimensions
x1 = np.random.normal(size=100)
x2 = np.random.normal(size=100)
x3 = x1 + x2
X = np.c_[x1, x2, x3]
from sklearn import decomposition
pca = decomposition.PCA()
pca.fit(X)
print(pca.explained_variance_)
# As we can see, only the 2 first components are useful
pca.n_components = 2
X_reduced = pca.fit_transform(X)
X_reduced.shape

独立成分分析：ICA

独立成分分析（ICA）选择合适的成分使得他们的分布载有最大的独立信息量。可以恢复非高斯独立信号：

# Generate sample data
time = np.linspace(0, 10, 2000)
s1 = np.sin(2 * time)  # Signal 1 : sinusoidal signal
s2 = np.sign(np.sin(3 * time))  # Signal 2 : square signal
S = np.c_[s1, s2]
S += 0.2 * np.random.normal(size=S.shape)  # Add noise
S /= S.std(axis=0)  # Standardize data
# Mix data
A = np.array([[1, 1], [0.5, 2]])  # Mixing matrix
X = np.dot(S, A.T)  # Generate observations
# Compute ICA
ica = decomposition.FastICA()
S_ = ica.fit_transform(X)  # Get the estimated sources
A_ = ica.mixing_.T
np.allclose(X,  np.dot(S_, A_) + ica.mean_)

五、联合起来

（1）管道（流水线）

我们已经知道了一些估测器（模型）能够转换数据，一些可以预测变量。我们也能够将其结合到一起：

from sklearn import linear_model, decomposition, datasets
from sklearn.pipeline import Pipeline
from sklearn.grid_search import GridSearchCV
logistic = linear_model.LogisticRegression()
pca = decomposition.PCA()
pipe = Pipeline(steps=[('pca', pca), ('logistic', logistic)])
digits = datasets.load_digits()
X_digits = digits.data
y_digits = digits.target
###############################################################################
# Plot the PCA spectrum
pca.fit(X_digits)
plt.figure(1, figsize=(4, 3))
plt.clf()
plt.axes([.2, .2, .7, .7])
plt.plot(pca.explained_variance_, linewidth=2)
plt.axis('tight')
plt.xlabel('n_components')
plt.ylabel('explained_variance_')
###############################################################################
# Prediction
n_components = [20, 40, 64]
Cs = np.logspace(-4, 4, 3)
#Parameters of pipelines can be set using ‘__’ separated parameter names:
estimator = GridSearchCV(pipe,
dict(pca__n_components=n_components,
logistic__C=Cs))
estimator.fit(X_digits, y_digits)
plt.axvline(estimator.best_estimator_.named_steps['pca'].n_components,
linestyle=':', label='n_components chosen')
plt.legend(prop=dict(size=12))

（2）使用特征联进行人脸识别

该实例使用的数据集是从“Labeled Faces in the Wild”节选预处理得到的。更为熟知的名字是LFW。

http://vis-www.cs.umass.edu/lfw/lfw-funneled.tgz（233 MB）

"""
===================================================
Faces recognition example using eigenfaces and SVMs
===================================================
The dataset used in this example is a preprocessed excerpt of the
"Labeled Faces in the Wild", aka LFW_: http://vis-www.cs.umass.edu/lfw/lfw-funneled.tgz (233MB)
.. _LFW: http://vis-www.cs.umass.edu/lfw/ Expected results for the top 5 most represented people in the dataset::
precision    recall  f1-score   support
Gerhard_Schroeder       0.91      0.75      0.82        28
Donald_Rumsfeld       0.84      0.82      0.83        33
Tony_Blair       0.65      0.82      0.73        34
Colin_Powell       0.78      0.88      0.83        58
George_W_Bush       0.93      0.86      0.90       129
avg / total       0.86      0.84      0.85       282
"""
from __future__ import print_function
from time import time
import logging
import matplotlib.pyplot as plt
from sklearn.cross_validation import train_test_split
from sklearn.datasets import fetch_lfw_people
from sklearn.grid_search import GridSearchCV
from sklearn.metrics import classification_report
from sklearn.metrics import confusion_matrix
from sklearn.decomposition import RandomizedPCA
from sklearn.svm import SVC
print(__doc__)
# Display progress logs on stdout
logging.basicConfig(level=logging.INFO, format='%(asctime)s %(message)s')
###############################################################################
# Download the data, if not already on disk and load it as numpy arrays
lfw_people = fetch_lfw_people(min_faces_per_person=70, resize=0.4)
# introspect the images arrays to find the shapes (for plotting)
n_samples, h, w = lfw_people.images.shape
# for machine learning we use the 2 data directly (as relative pixel
# positions info is ignored by this model)
X = lfw_people.data
n_features = X.shape[1]
# the label to predict is the id of the person
y = lfw_people.target
target_names = lfw_people.target_names
n_classes = target_names.shape[0]
print("Total dataset size:")
print("n_samples: %d" % n_samples)
print("n_features: %d" % n_features)
print("n_classes: %d" % n_classes)
###############################################################################
# Split into a training set and a test set using a stratified k fold
# split into a training and testing set
X_train, X_test, y_train, y_test = train_test_split(
X, y, test_size=0.25)
###############################################################################
# Compute a PCA (eigenfaces) on the face dataset (treated as unlabeled
# dataset): unsupervised feature extraction / dimensionality reduction
n_components = 150
print("Extracting the top %d eigenfaces from %d faces"
% (n_components, X_train.shape[0]))
t0 = time()
pca = RandomizedPCA(n_components=n_components, whiten=True).fit(X_train)
print("done in %0.3fs" % (time() - t0))
eigenfaces = pca.components_.reshape((n_components, h, w))
print("Projecting the input data on the eigenfaces orthonormal basis")
t0 = time()
X_train_pca = pca.transform(X_train)
X_test_pca = pca.transform(X_test)
print("done in %0.3fs" % (time() - t0))
###############################################################################
# Train a SVM classification model
print("Fitting the classifier to the training set")
t0 = time()
param_grid = {'C': [1e3, 5e3, 1e4, 5e4, 1e5],
'gamma': [0.0001, 0.0005, 0.001, 0.005, 0.01, 0.1], }
clf = GridSearchCV(SVC(kernel='rbf', class_weight='auto'), param_grid)
clf = clf.fit(X_train_pca, y_train)
print("done in %0.3fs" % (time() - t0))
print("Best estimator found by grid search:")
print(clf.best_estimator_)
###############################################################################
# Quantitative evaluation of the model quality on the test set
print("Predicting people's names on the test set")
t0 = time()
y_pred = clf.predict(X_test_pca)
print("done in %0.3fs" % (time() - t0))
print(classification_report(y_test, y_pred, target_names=target_names))
print(confusion_matrix(y_test, y_pred, labels=range(n_classes)))
###############################################################################
# Qualitative evaluation of the predictions using matplotlib
def plot_gallery(images, titles, h, w, n_row=3, n_col=4):
"""Helper function to plot a gallery of portraits"""
plt.figure(figsize=(1.8 * n_col, 2.4 * n_row))
plt.subplots_adjust(bottom=0, left=.01, right=.99, top=.90, hspace=.35)
for i in range(n_row * n_col):
plt.subplot(n_row, n_col, i + 1)
plt.imshow(images[i].reshape((h, w)), cmap=plt.cm.gray)
plt.title(titles[i], size=12)
plt.xticks(())
plt.yticks(())
# plot the result of the prediction on a portion of the test set
def title(y_pred, y_test, target_names, i):
pred_name = target_names[y_pred[i]].rsplit(' ', 1)[-1]
true_name = target_names[y_test[i]].rsplit(' ', 1)[-1]
return 'predicted: %s\ntrue:      %s' % (pred_name, true_name)
prediction_titles = [title(y_pred, y_test, target_names, i)
for i in range(y_pred.shape[0])]
plot_gallery(X_test, prediction_titles, h, w)
# plot the gallery of the most significative eigenfaces
eigenface_titles = ["eigenface %d" % i for i in range(eigenfaces.shape[0])]
plot_gallery(eigenfaces, eigenface_titles, h, w)
plt.show()

预测：

特征脸：

数据集中最有代表性的五个人脸的期望结果：

precision recall f1-score support

Gerhard_Schroeder 0.91 0.75 0.82 28

Donald_Rumsfeld 0.84 0.82 0.83 33

Tony_Blair 0.65 0.82 0.73 34

Colin_Powell 0.78 0.88 0.83 58

George_W_Bush 0.93 0.86 0.90 129

avg / total       0.86      0.84      0.85       282

（3）开放性问题：股票市场结构

我们是否可以根据给定的时间帧预测股票的价格变化。
[学习一个图结构]

六、寻求帮助

###（1）项目邮件列表

如果你碰到scikit-learn的BUG或者文档中需要澄清声明的部分，请放心大胆的在邮件列表里询问[maillist]

(2）问答（Q&A）机器学习从业者参与的社区

Metaoptimize/QA:

一个 机器学习、自然语言处理和其他数据分析方面讨论的论坛（类似针对开发者的Stackoverflow）:http://metaoptimize.com/qa

一个比较容易开始参与的讨论：good freely available textbooks on machine learning（机器学习方面优秀的免费电子书）

Quora.com:

Quora 有一个关于机器学习相关的问题主题，也有很多有趣的讨论：http://quora.com/Machine-learning

浏览一下最佳问题的部分，例如：What are some good resources for learning about machine learning(关于机器学习的优秀资源有哪些)

---斯坦福的 Andrew Ng教授 教授的 关于机器学习的优秀在线免费课程

{网易公开课有，搜一下机器学习就可以了}

---一个更倾向于人工智能（AI)的优秀在线课程:
http://www.udacity.com/overview/Course/cs271/CourseRev/1