一、基础知识：
1. 数据预处理：
实际应用时经常需要处理大量原始数据，为了让相应软件可以使用这些数据，往往需要对数据预处理。
1）均值移除Mean removal：
通过会把每个特征的平均值移除，以保证均值为0，这样可以消除特征彼此之间的偏差。
from sklearn import preprocessing
datastand=preprocessing.scale(data)
print("\nMean=",datastand.mean(axis=0))
print("Std deviation=",datastand.std(axis=0))
数据处理后特征均值几乎为0，标准差则为1。
2）范围缩放Scaling：
数据点中每个特征的数值范围可能变化很大，因此有时需要将特征的数值范围缩放到合理的大小。
datascaler=preprocessing.MinMaxScaler(feature_range=(0,1))
datascaled=datascaler.fit_transform(data)
print("\nMin max scaled data=",datascaled)
范围缩放后，所有数据点的特征数值都位于指定的数值范围内。
3）归一化Normalization：
数据归一化用于对特征向量的值进行调整，以保证每个特征向量的值都缩放到相同的数值范围。机器学习中最常用的归一化形式是将特征向量调整为L1范数，使特征向量的数值之和为1。
datanormal=preprocessing.normalize(data,norm='l1')
print('\nL1 normalized data=',datanormal)
这个方法经常用于确保特征点数据处于同一数量级，提高不同特征数据的可比性。
4）二值化Binarization：
二值化用于将数值特征向量转换为布尔型向量。
databinarized=preprocessing.Binarizer(threshold=1.4).transform(data)
print('\nBinarized data=',databinarized)
5）独热编码One-Hot Encoding：
通常，需要处理的数据都是稀疏的，散乱地分布在空间。独热编码可看作是一种收紧特征向量的工具，它把特征向量的每个特征与特征的非重复总数相对应，这样可以更加有效地表示数据空间。比如，编码器遍历每个特征向量的第n个特征，然后进行非重复计数，如果非重复计数的值是K，那么就把这个特征转换为只有一个值是1其他值都是0的K维向量。
encoder=preprocessing.OneHotEncoder()
encoder.fot([0,2,1,12],[1,3,5,3],[2,3,2,12],[1,2,4,3]])
encodedvector=encoder.transform([[2,3,5,3]]).toarray()
print('\nEncoded vector=',encodedvector)
2. 标记编码方法：
在监督学习中，经常需要处理各种各样的标记。这些标记可能是数字，有可能是单词。如果标记是数字，算法可以直接使用，但很多情况下标记需要以人可理解的形式存在，因此通常会用单词标记训练集。标记编码就是要把单词标记转换成数值形式，让算法懂得如何操作标记。
⑴定义编码器label encoder：
label_encoder=preprocessing.LabelEncoder()
⑵创建标记：
inputclass=['audi','ford','audi','toyota','ford','bmw']
⑶为标记编码：
label_encoder.fit(inputclass)
print('\nClass mapping:')
for i,item in enumerate(label_encoder.classer_):
print(item,'-->',i)
这样，单词被转换成从0开始的索引值。然后，如果遇到一组标记，就可以进行如下处理：
labels=['toyota','for','audi']
labelencoded=label_encoder.transform(labels)
print('\nLabels=',labels)
print('Encoded labels=',list(labelencoded)
这种方法比手动进行单词与数字的编码要简单，还可以通过数字反转会单词的功能检查结果：
labelencoded=[2,1,0,3,1]
labeldecoded=label_encoder.inverse_transform(labelencoded)
print('\nEncoded labels=',labelencoded)
print('Decoded labels=',list(labeldecoded)
3. 构建对数据集中每个数据处理的函数：
1）将数组的每一个元素加3：
def add3(input_array):
return map(lambda x:x+3, input_array)
2）将数组的每一个元素乘以2：
def mul2(input_array):
return map(lambda x:x*2, input_array)
3）将数组的每一个元素减去5：
def sub5(input_array):
return map(lambda x:x-5, input_array)
4）定义函数组合器：
def function_composer(*args):
return reduce(lambda f,g:lambda x:f(g(x)),args)
5）单行代码的函数组合器：
func_compised=function_composer(sub5, nul2, add3)
二、回归器：
1. 线性回归器：
回归是估计输入数据与连续值输出数据之间关系的过程。数据通常是实数形式，目标是估计满足输入到输出映射关系的基本函数。线性回归用输入变量的线性组合来估计基本函数。
线性回归的目标是提取的输入变量与输出变量的关联模型，这就要求实际输出与线性方程预测的输出的残差平方和最小化，这种方法被称为普通最小二乘法OLS（Ordinary Least Squares）。
线性回归的主要优点是方程简单，而使用非线性回归的曲线方程拟合速度会慢很多。
1）创建线性回归器：
如果已经创建了数据文件data.txt，文件中用逗号分隔符分隔字段，第1个字段是输入值，第2个字段是与输入值相对应的输出值。
filename=sys.argv[1]
X=[];y=[]
with open(filename,'r') as f:
for line in f.readlines():
xt,yt=[floar(i) for i in line.split(',')]:
X.append(xt);y.append(yt)
代码中把输入数据加载到变量X和y，其中X是数据，y是标记。在for循环中，解析每行数据，用逗号分隔，然后将字段转化为浮点数，并分别保存到变量X和y中。
建立机器学习模型时，把数据分成两组，训练数据集training dataset与测试数据集resting dataset。训练数据集用来建立模型，测试数据集用来验证模型对未知数据的学习效果。
num_training=int(0.8*len(X))
num_test=len(X)-num_training
X_train=np.array(X[:num_training]).reshape((num_training,1))
y_train=np.array(y[:num_training])
X_test=np.array(X[num_training:]).reshape((num_test,1))
y_test=np.array(y[num_training:])
这里用80%的数据作为训练数据集，其余20%的数据作为测试数据集。
接下来创建一个回归器对象：
from sklearn import linear_model
linear_regressor=linear_model.LinearRegression()
linear_regressor.fit(X_train,y_train)
然后利用训练数据集训练线性回归器，向fit方法提供输入数据即可训练模型。可以通过绘图来直观看到结果。
import matplotlib.pyplot as plt
y_train_pred=linear_regressor.predict(X_train)
plt.figure()
plt.scatter(X_train, y_train, color='green')
plt.plot(X_train, y_train_pred, color='black', linewidth=4)
plt.title('Training data')
plt.show()
接下来要用模型对测试数据集进行预测，并画出图形：
y_test_pred=linear_regressor.predice(X_test)
plt.scatter(X_test, y_test, color='green')
plt.plot(X_test, y_test_pred, color='black', linewidth=4)
plt.title('Test data')
plt.show()
2）计算回归准确性：
建立的回归器，然后要评价拟合效果，用误差表示实际值与模型预测值之间的差值。回归器可以使用不同的指标进行衡量，部分指标有：
·平均绝对误差mean absolute error：是给定数据集的所有数据点的绝对误差平均值
·均方误差mean squared error：给定数据集的所有数据点的误差的平方的平均值
·中位数绝对误差median absolute error：给定数据集的所有数据点的误差的中位数
·解释方差分explained variance score：用于衡量模型对数据集波动的解释能力的分数
·R方得分R2 score：是指确定性相关系数，用于衡量模型对未知样本预测的效果
scikit-learn中的一个模块提供了计算所有指标的功能：
import sklearn.metrics as sm
print('Mean absolute error=', round(sm.mean_absolute_error(y_test, y_test_pred), 2)
print('Mean aquared error=', round(sm.mean_squared_error(y_test, y_test_pred), 2)
print('Median absolute error=', round(sm.median_absolute_error(y_test, y_test_pred), 2)
print('Explained variance score=',round(sm.explained_variance_score(y_test, y_test_pred), 2)
print('R2 score=',round(sm.re_score(y_test, y_test_pred), 2)
一般只选择一两个指标来评估模型，通常是尽量保证均方误差最低，且方差分最高。
3）保存模型数据：
模型训练结束后，把模型保存为文件，下次再使用时只需要简单地加载就可以。
import cPickle as pickle
output_model_file='saved_model.pkl'
with open(output_model_file, 'w') as f:
pickle.dump(linear_regressor, f)
回归模型会保存在saved_model.plk文件中。以后可以加载并使用：
with open(output_model_file, 'r') as f:
model_linregr=pickle.load(f)
y_test_pred_new=model_linregr.predict(X_text)
print('\nNew mean absolute error=', round(sm.mean_absolute_error(y_test, y_test_pred_new), 2)
代码中将回归模型从Pickle文件加载到model_linregr变量中。
2. 岭回归器：
线性回归对异常值比较敏感，实际收集的数据中经常会遇到错误的结果，误差值绝对值很大就能破坏整个模型。为了避免此问题，引入正则化项的系数作为阈值来消除异常值的影响，这种方法称为岭回归。
1）创建岭回归器：
从data_nulti.txt文件中加载数据，其中每一行都包含多个值，除最后一个值外，前面所有数值构成输入特征向量。使用参数初始化岭回归器：
ridge_regressor=linear_model.Ridge(alpha=0.01, fit_intercept=True, max_iter=10000)
其中alpha参数控制回归器的复杂程序，当alpha趋于0时就是普通线性回归器。如果希望模型对异常值不敏感，就需要设置一个较大的值。
2）训练岭回归器：
ridge_regressor.fit(X_train, y_train)
y_test_pred_ridge=ridge_regressor.predict(X_test)
print('Mean absolute error=', round(sm.mean_absolute_error(y_test, y_test_pred_ridge), 2)
print('Mean aquared error=', round(sm.mean_squared_error(y_test, y_test_pred_ridge), 2)
print('Median absolute error=', round(sm.median_absolute_error(y_test, y_test_pred_ridge), 2)
print('Explained variance score=',round(sm.explained_variance_score(y_test, y_test_pred_ridge), 2)
print('R2 score=',round(sm.re_score(y_test, y_test_pred_ridge), 2)
运行代码检查误差指标。
3. 多项式回归器：
线性回归模型只能把输入数据拟合成直线，而多项式回归模型通过拟合多项式非常可以获得曲线，从而提高模型的准确性。
import sklearn.preprocessing import PolynomialFeatures
polynomial=PolynomialFeatures(degree=3)
代码将曲线的多项式的次数初始值设为3.然后用数据点来计算多项式的参数：
X_train_transformed=polynomial.fit_transform(X_train)
其中，X_train_transformed表示多项式形式的输入。接下来用文件中的第一个数据点来检查多项式模型是否能准确预测：
datapoint=[0.39, 2.78, 7.11]
poly_datapoint=polynomial.fit_transform(datapoint)
poly_linear_model=linear_model.LinearRegression()
poly_linear_model.fit(X_train_transformed, y_train)
print('\nLinear regression:', linear_regressor.predict(datapoint)[0]
print('\nPolynimial regression:', poly_linear_model.predict(poly_datapoint)[0]
多项式回归模型计算变量数据点的值恰好就是输入数据文件中的第一行数据值。再用线性回归模型测试以下，比较结果。如果想要数据更接近实际输出值，就需要增加多项式的次数。可以将多项式的次数加到10：
polynomial=PolynomialFeatures(degree=10)
增加多项式模型的复杂性拟合速度会变慢。
4. 决策树回归器：
决策树是一个树状模型，每个节点都做出一个决策，从而影响最终结果。叶子节点表示输出数值，分支表示根据输入特征做出的中间决策。AdaBoost算法是指自适应增强adaptive boosting算法，将不同版本的算法结果进行组合，用加权汇总的方式得到最终结果，这称为弱学习器weak learners。AdaBoost算法在每个阶段获取的信息都会反馈到模型中，这样学习器就可以在后一阶段重点训练难以分类的样本，可以增强系统的准确性。
首先使用AdaBoost算法对数据集合进行回归拟合，再计算误差，然后根据误差评估结果，用同样的数据集重新拟合。可以把这些看作是回归器的调优过程，直到达到预期的准确性。
import numpy as np
from sklearn.tree import DecisionTreeRegressor
from sklearn.ensemble import AdaBoostRegressor
from sklearn import datasets
from sklearn.metrics import mean_squared_error, explained_variance_score
from sklearn.utils import shuffle
import matplotlib.pyplot as plt
输入数据。这是一个房价数据，每个数据点由13个影响房价的参数构成。
hdata=datasets.load_boston()
接下来把数据与输出结果分成不同的变量，通过shuffle函数把数据的顺序打乱：
X,y=shuffle(hdata.data, hdata.target, random_state=7)
其中random_state用来控制打乱数据方法。然后把数据分成训练数据集和测试数据集，其中80%用于训练，剩余20%用于测试。
num_training=int(0.8*len(X))
X_train, y_train=X[:num_training], y[:num_training]
X_test, y_test=X[num_training:], y[num_training:]
这样就可以拟合决策树回归模型了，选一个最大深度为4的决策树：
dt_regressor=DecisionTreeRegressor(man_depth=4)
dt_regressor.fit(X_train, y_train)
再用带AdaBoost算法的决策树回归模型进行拟合：
ab_regressor=AdaBoostRegressor(DecisionTreeRegressor(max_depth=4), n_estimators=400, random_state=7)
ab_regressor.fit(X_train, y_train)
接下来评价决策树回归器的训练效果：
y_pred_dt=dt_regressor.predict(X_test)
mse=mean_squared_error(y_test, y_pred_dt)
evs=explained_variance_score(y_test, y_pred_dt)
print('\nMean squared error=', round(mse, 2))
print('Explained variance score=', round(evs, 2))
再评价AdaBoost算法改善：
y_pred_ab=ab_regressor.predict(X_test)
mse=mean_squared_error(y_test, y_pred_ab)
evs=explained_variance_score(y_test, y_pred_ab)
print('\nMean squared error=', round(mse, 2))
print('Explained variance score=', round(evs, 2))
5. 计算特征的相对重要性：
对一些问题，会使用多个特征，它们对模型都有贡献，但对结果的贡献并不一样。
plot_feature_importances(dt_regressor.feature_importances_, 'Decision Tree regressor', hdata.feature_names)
plot_feature_importances(ab_regressor.feature_importances_, 'AdaBoost regressor', hdata.feature_names)
回归对象有一个feature_importances_方法会得到每个特征的相对重要性。然后就可以画图：
def plot_feature_importances(feature_importances, title, feature_names):
feature_importances=100.0*(feature_importances/max(feature_importances))
index_sorted=np.flipud(np.argsort(feature_importances))
pos=np.arange(index_sorted.shape[0])+0.5
plt.figure()
plt.bar(pos, feature_importances[index_sorted], align='center')
plt.xticks(pos, feature_names[index_sorted])
plt.ylabel('Relative Importance')
plt.title(title)
plt.show()
代码中从feature_importances_方法中取值，然后把数值放大到0~100的范围内。
6. 随机森林回归器：
随机森林回归器random forest regressor是一个决策树集合，就是用一组由数据集的若干子集构建的决策树构成，再用决策树平均值改善整体学习效果。
使用bike.csv数据文件，共16列，前两列是序列号与日期，分析时不用；最后3列数据是不同类型的输出结果；最后1列是14与15列的和，因此也不用。
import csv
from sklearn.ensemble import RandomForestRegressor
from housing import plot_feature_importances
需要处理CSV文件，所以加入csv程序包来读取，还需要自定义数据集加载函数：
def load_dataset(filename):
file_reader=csv.reader(open(filename, 'rb'), delimiter=',')
X,y=[],[]
for row in file_reader:
X.append(row[2:13])
y.append(row[-1])
feature_names=np.array(X[0]) # 提取特征名称
return np.array(X[1:]).astype(np.float32), np.array(y[1:]).astype(np.float32), feature_names
此函数中从CSV读取所有数据，把特征名称从其中分离出来并作为函数返回值。然后要打乱数据顺序。
X,y=feature_names=load_dataset(sys.argv[1])
X,y=shuffle(X, y, random_state=7)
需要将数据集分成训练集和测试集，以90%数据用于训练：
num_training=int(0.9*len(X))
X_train, y_train=X[:num_training], y[:num_training]
X_test, y_test=X[num_training:], y[num_training:]
然后训练回归器：
rf_regressor=RandomForestRegressor(n_estimators=1000, max_depth=10, min_samples_split=1)
rf_regressor.fit(X_train, y_train)
其中，n_estimators是指评估器的数量，表示随机森林需要使用的决策数数量；max_depth是指每个决策树的最大深度；min_samples_split是指决策树分裂一个节点需要用到的最小数据样本量。
评价随机森林回归器的训练效果：
y_pred=rf_regressor.predict(X_test)
mse=mean_squared_error(y_test, y_pred)
evs=explained_variance_score(y_test, y_pred)
print('\nMean squared error=', round(mse, 2))
print('Explained variance score=', round(evs, 2))
已经有了画出特征重要性条形图的函数plot_feature_importances，接下来直接调用：
plot_feature_importances(rf_regressor.feature_importance_, 'Random Forst regressor', feature_names)
三、分类器classifier：
在机器学习领域，分类是指利用数据的特性将其分为若干类型。分类器可以是实现分类功能的任意算法，最简单的分类器就是简单的数学函数。分类器有二元binary分类器和多元multiclass分类器，解决分类问题的数据手段都倾向于解决二元分类问题，可以通过扩展而解决多元分类问题。
1. 简单分类器：
数据样本：X=np.array([[3,1],[2,5],[1,8],[6,4],[5,2],[3,5],[4,7],[4,-1]])
为这些数据点分配一些标记：y=[0,1,1,0,0,1,1,0]
因为只有两类，y列表中只包含0和1.一般情况下，如果有N个类，y取值范围就是从0到N-1。然后按照类型标记把样本数据分成两类：
class_0=np.array([X[i] for i in tange(len(X)) if y[i]==0])
class_1=np.array([X[i] for i in tange(len(X)) if y[i]==1])
使用图像方式画出：
plt.figure()
plt.scatter(class_0[:0], class_0[:1], color='black', marker='s')
plt.scatter(class_1[:0], class_1[:1], color='black', marker='x')
这是散点图scatterplot，用方块和叉表示两类数据。
在两类数据间画一条分割线，只要用直线方程就可以：
line_x=range(10)
line_y=line_x
在图像中显示：
plt.plot(line_x, line_y, color='black', linewidth=3)
plt.show()
上述示例的分类器，如果输入点(a,b)的a大于等于b就属于类型class_0；否则就属于class_1。这是一个线性分类器，分割线是一条直线。如果分割线是一条曲线，就像非线性分类器。
2. 逻辑回归分类器：
给定一组数据点，建立一个可以在类之间绘制线性边界的模型。逻辑回归可以对训练数据派生的一组方程进行求解来提取边界。假设有9组数据，分为3类。
X=mp.array([[4,7],[3.5,8],[3.1,6.2],[0.5,1],[1,2],[1.2,1.9],[6,2],[5.7,1.5],[5.4,2.2]])
y=np.array([0,0,0,1,1,1,2,2,2])
初始化逻辑回归分类器：
classifier=linear_model.LogisticRegression(solver='liblinear',C=100)
函数中有多个参数，其中最重要的是solver和C，solver用于设置求解系统方程的算法模型，C表示正则化强度，数值越小表示正则化强度越强。
训练分类器：
classifier.fit(X,y)
画出数据点和边界：
plot_classifier(classifier,X,y)
其中使用了自定义画图函数：
def plot_classifier(classifier,X,y):
x_min,x_max=min(X[:,0])-1.0,max(X[:,0)+1.0
y_min,y_max=min(X[:,1])-1.0,max(X[:,1)+1.0
在图形中要使用的数值范围通常是从最小值到最大值，并增加了一个余量。
为了画出边界，需要利用一组网格grid数据求出方程的值，然后把边界画出来：
step_size=0.01
x_values, y_values=np.meshgrid (np.arange (x_min, x_max, step_size), np.arange (y_min, y_max, step_size))
变量x_values和y_values包含了求解方程数值的网格点。
计算出分类器对所有数据点的分类结果：
mesh_output=classifier.predict(np.c_[x_values.ravel(), y_values.ravel()])
mesh_output=mesh_output.reshape(x_values.shape)
用彩色区域画出各类型的边界：
plt.figure()
plt.pcolormesh(x_values, y_values, mesh_output, cmap=plt.cm.gray)
这是一个三维画图器，也可以画二维数据点，可以用色彩清单表示不同区域的相关属性。
接着把训练数据画在图上：
plt.scatter(X[:,0], X[:,1], c=y, s=80, edgecolors='black', linewidth-1, cmap=plt.cm.Paired)
plt.xlim(x_values.min(), x_values.max())
plt.ylim(y_values.min(), y_values.max())
plt.xticks((np.arange(int(min(X[:,0])-1), int(max(X[:,0])+1), 1.0)))
plt.yticks((np.arange(int(min(X[:,0])-1), int(max(X[:,0])+1), 1.0)))
plt.show()
其中，plt.scatter把数据点画在二维图上，X[:,0]表示0轴坐标值，X[:,1]表示1轴坐标值，c=y表示颜色的使用顺序，cmap表示映射的颜色表。坐标值的取值范围由plt.xlim和plt.ylim确定，标记坐标轴的数值，需要使用plt.xticks和plt.yticks，希望刻度是整数，所以使用int()取整。
可以通过调整参数C看对模型的影响，C表示对分类错误的惩罚值，设为1.0、10000会得到不同的结果。随着参数C的不断增大，分类错误的惩罚值越高，各个类型的边界更优。
3. 朴素贝叶斯分类器：
朴素贝叶斯分类器是用贝叶斯定理进行建模的监督学习分类器。首先要导入程序包：
from sklearn.naive_bayes import GaussianNB
from logistic_regression import plot_classifier
数据来自txt文件，每一行数据都是由逗号分隔的数值。
imput_file='data.txt'
X=[];y=[]
with open(input_file, 'r') as f:
for line in f.readlines():
data=[float(x) for x in line.split(',')]
X.append(data[:-1]);y.append(data[-1])
X=np.array(X);y=np.array(y)
上面代码已经把数据和标记分别加载到变量X和y中了。需要将数据分为训练集和测试集：
X_train, X_test, y_train, y_test = cross_validation.train_test_split (X, y, test_size=0.25 , random_state=5)
然后建立朴素贝叶斯分类器：
classifier_gaussiannb=GaussianNB()
classifier_gaussiannb.fit(X_train,y_train)
y_pred=classifier_gaussiannb.predict(X_train, y_train)
GaussianNB()函数指定了正态分别贝叶斯模型Gaussian Naive Bayes model。然后计算分类器的准确性：
accuracy=100.0*(y_test==y_test_pred).sum()/X_test.shape[0]
print('Accuracy of the classifier=', round(accuracy, 2), '%')
接着画出数据点和边界：
plot_classifier(classifier_gaussiannb, x_test, y_test)
4. 机器学习的指标：
处理机器学习模型时，通常关心3个指标，精度precision、召回率recall和F1得分F1 score。精度指被分类器正确分类的样本数量占分类器总分类样本数量的百分比；召回率是指被正确分类的样本数量占某分类总样本数量的百分比。
精度和召回率是背反的，需要保持两个指标能同时处于合理值，引入了F1得分指标，这是精度和召回率的合成指标，是调和均值：
F1得分=2x精度x召回率/(精度+召回率)
计算上面示例的指标：
num_validations=5
accuracy = cross_validation.cross_val_score (classifier_gaussiannb, X, y, scoring='accuracy', cv=num_validations)
print('Accuracy: '+str(round(100*accuracy.mean(), 2))+'%'
f1 = cross_validation.cross_val_score (classifier_gaussiannb, X, y, scoring='f1_weighted', cv=num_validations)
print('F1: '+str(round(100*f1.mean(), 2))+'%'
precision = cross_validation.cross_val_score (classifier_gaussiannb, X, y, scoring='precision_weighted', cv=num_validations)
print('Precision: '+str(round(100*precision.mean(), 2))+'%'
recall = cross_validation.cross_val_score (classifier_gaussiannb, X, y, scoring='recall_weighted', cv=num_validations)
print('Recall: '+str(round(100*recall.mean(), 2))+'%'
5. 混淆矩阵可视化：
混淆矩阵confusion matrix是理解分类模型性能的数据表。
| Predicted class0 | Predicted class1 | Predicted class2 | |
| True class0 | 45 | 4 | 3 |
| True class1 | 11 | 56 | 2 |
| True class2 | 5 | 6 | 49 |
图中可以看出不同类型的分类数据。理想情况下，希望矩阵非对角线元素都是0，这是最完美的分类结果。
from sklearn.metrics import confusion_matrix
y_true=[1,0,0,2,1,0,3,3,3]
y_pred=[1,1,0,2,1,0,1,3,3]
confusion_mat=confusion_matrix(y_true, y_pred)
plot_confusion_matrix(confusion_mat)
其中使用了自定义的混淆矩阵画图函数：
def plot_confusion_matrix(confusion_mat)
plt.imshow(confusion_mat, interpolation='nearest', cmap=plt.cm.Paired)
plt.title('Confusion matrix')
plt.colorbar()
tick_marks=np.arange(4)
plt.xticks(tick_marks, tick_marks)
plt.yticks(tick_marks, tick_marks)
plt.ylabel('True label')
plt.xlabel('Predicted label')
plt.show()
其中使用imshow函数画混淆矩阵，参数tick_marks的取值范围0~3，因为数据集中有4个标记类型。画出的图中，对角线颜色很亮，而黑色区域表示0，灰色区域表示分类错误的样本量。
6. 提取性能报告：
from sklearn.metrics import classification_report
y_true=[1,0,0,2,1,0,3,3,3]
y_pred=[1,1,0,2,1,0,1,3,3]
target_names=['Class-0','Class-1','Class-2','Class-3']
print(classification_report(y_true, y_pred, target_names=target_names))
7. 生成学习曲线：
学习曲线可以帮助理解训练数据集大小对机器学习模型的影响。
from sklearn.learning_curve import learning_curve
classifier=RandomForestClassifier(random_state=7)
parameter_grid=np.array([200,500,800,1100])
train_sizes, train_scores, validation_scores = learning_curve (classifier, X, y, train_sizes =parametre_grid, cv=5)
print('nTraining scores: \n', train_scores)
print('\nValidation scores:\n', validation_scores)
代码中分别用200、500、800、1100的训练数据集的大小测试模型的性能指标，将cv参数设置为5，就是用五折交叉验证。
然后把数据画成图形：
plt.figure()
plt.plot(parameter_grid, 100*np.average(train_scores, axis=1), color='black')
plt.title('Learning curve')
plt.xlable('Number of training samples')
plt.ylabel('Accuracy')
plt.show()
四、预测建模：
预测建模Predictive modeling是一种用来预测系统未来行为的分析技术，它由一群能够识别独立输入变量与反馈目标关联关系的算法构成。根据预测值创建一个数学模型，然后用这个模型去预测未来的趋势。
预测建模中，需要收集已知的响应数据来训练模型。一旦模型建成，就可以用一些指标来检验，然后用它预测未来值。
1. 用SVM建立线性分类器：
SVM(Support Vector Machine)是用来构建分类器和回归器的监督学习模型。SVM通过数学方程组求解，可以找出两组数据之间的最佳分割边界。分割数据为训练集和测试集：
from sklearn import cross_validation
from sklearn.svm import SVC
X_train, x_test, y_train, y_test = cross_validation.train_test_split (X, y, test_size=0.25, random_state=5)
用线性核函数初始化一个SVM对象：
params={'kernel':'linear'}
classifier=SVC(**params)
训练SVM分类器：
classifier.fit(X_train, y_train)
可视化：
utilities.plot_classifier(classifier, X_train, y_train, 'Training dataset')
plt.show()
测试分类器：
y_test_pred=classifier.predict(X_test)
utilities.plot_classifier(classifier, X_test, y_test, 'Test dataset')
plt.show()
接着计算训练集的准确性：
from sklearn.metrics import classification_report
target_names=['Class-'+str(int(i)) for i in set(y)]
print(classification_repart(y_train, classifier.predict(X_train), target_names=target_names))
分类器为测试集生成的分类报告：
print(classifier_report(y_test, y_test_pred, target_names=target_names))
2. 用SVM建立非线性分类器：
SVM为建立非线性分类器提供了许多选项，需要使用不同的核函数建立。常用的是多项式函数和径向基函数RBF（Radial Basic Function）。
用一个多项式核函数建立非线性分类器：
params={'kernel':'poly','degree':3}
这里用了三次多项式方程。可以增加方程次数，但训练时间更长。
也可以使用径向基函数建立非线性分类器：
params={'kernel':'rbf'}
有时，某一类型的数据点数量会比其他类型的多很多，这时训练的分类器就会有偏差。
params={'kernel':'linear','class_weight':'auto'}
参数class_weight的作用是统计不同类型数据点的数量，调整权重，让类型不平衡问题不影响分类效果。
3. 提取置信度：
当一个新的数据点被分类为某个已知类别时，可以训练SVM来计算出输出类似的置信度。
input_datapoints=np.array([[2,1.5],[8,9],[4.8,5.2],[4,4],[2.5,7],[7.6,2],[5.4,5.9]])
测量数据点到边界的距离：
print('\nDistance from the boundary:')
for i in input_datapoints:
print(i, '-->', classifier.decision_function(i)[0])
概率输出是一种将不同类型的距离度量转换成概率度量的方法。测量置信度：
params={'kernel':'rbf','probability':True}
classifier=SVC(**params)
参数probability用于指示SVM训练时要计算出概率。然后训练分类器：
classifier.fit(X_train, y_train)
计算输入数据点的置信度：
print('\nConfidence measure:')
for i in input_datapoints:
print(i, '-->', classifier.predict_proba(i)[0])
查看数据点与边界的位置：
utilities.plot_classifier(classifier, input_datapoints, [0]*len(input_datapoints), 'Input datapoints', 'True')
4. 寻找最优超参数：
parameter_grid=[{'kernel':['linear'], 'C':[1,10,50,600]},
{'kernel':['poly'], 'degree':[2,3]},
{'kernel':['rbf'], 'gamma':[0.01,0.001], 'C':[1,10,50,600]},]
定义需要使用的指标：
metrics=['precision', 'recall_weighted']
为每个指标搜索最优超参数：
for metric in metrics:
classifier=grid_search.GridSearchCV(svm.SVC(C=1), parameter_grid, cv=5, scoring=metric)
classifier.fit(X_train, y_train)
查看指标得分：
for params, avg_score, _ in classifier.grid_scores_:
print(params, '-->', round(avg_score, 3)
输出最好的参数集，模型将会尝试各种不同参数的组合来搜索最佳参数：
print('\nHighest scoring parameter set:', classifier.best_params_)
测试集上测试：
y_true, y_pred=y_test, classifier.predict(X_test)
print(classification_report(y_true, y_pred)
5. 建立事件预测器：
1）示例：预测大楼进出人数
数据每一行都由6个逗号分割的字符串组成，预测大楼举办的活动。
X=[];count=0
with open(input_file,'r') as f:
for line in f.readlines():
data=line[:-1].split(',')
X.append([data[0]+data[2:])
X=np.array(X)
label_encoder=[] # 字符串格式转换成数值格式
X_encoded=np.empty(X.shape)
for i,item in enumerate(X[0]):
if item.isdigit():
X_encoded[:i]=X[:,i]
else:
label_encoder.append(preprocessing.LabelEncoder())
X_encoded[:i]=label_encoded[-1].fit_transform(X[:i])
X=X_encoded[:,:-1].astype(int)
y=X_encoded[:,-1].astype(int)
params={'kernel':'rbf', 'probability':True, 'class_weight':'auto'} # 建立SVM模型
classofier=SVC(**params)
classifier.fit(X,y)
from sklearn import cross_validation # 交叉验证
accuracy=cross_validation.cross_val_score(classifier, X, y, scoring='accuracy', cv=3)
print('Accuracy of the classifier: '+str(round(100*accuracy.mean(), 2))+'%')
input_data=['Tuesday', '12:30:00', '21', '23'] # 对单一数据示例进行编码测试
input_data_encoded=[-1]*len(input_data]
count=0
for i,item in enumerate(input_data):
if item.isdigit():
input_data_encoded[i]=int(input_data[i])
else:
input_data_encoded[i]=int(label_encoder[count].transform(input_data[i]))
count=count+1
imput_data_encoded=np.array(input_data_encoded)
output_class=classifier.predict(input_data_encoded) # 为特定数据点预测并打印输出结果
print('Output class: ', label_encoded[-1].inverse_transform(output_class)[0]
2）估算交通流量：
X=[];count=0
with open(input_file,'r') as f:
for line in f.readlines():
data=line[:-1].split(',')
X.append(data)
X=np.array(X)
label_encoder=[] # 字符串格式转换成数值格式
X_encoded=np.empty(X.shape)
for i,item in enumerate(X[0]):
if item.isdigit():
X_encoded[:i]=X[:,i]
else:
label_encoder.append(preprocessing.LabelEncoder())
X_encoded[:i]=label_encoded[-1].fit_transform(X[:i])
X=X_encoded[:,:-1].astype(int)
y=X_encoded[:,-1].astype(int)
params={'kernel':'rbf', 'C':10.0, 'epsilon':0.2} # 建立SVR模型
regressor=SVR(**params)
repressor.fit(X,y)
omport sklearn.metrics as sm # 交叉验证
y_pred=regressor.predict(X)
print('Mean absolute error= '+round(sm.mean_absolute_error(y, y_pred),2))
input_data=['Tuesday', '13:35', 'San Francisco', 'yes'] # 对单一数据示例进行编码测试
input_data_encoded=[-1]*len(input_data]
count=0
for i,item in enumerate(input_data):
if item.isdigit():
input_data_encoded[i]=int(input_data[i])
else:
input_data_encoded[i]=int(label_encoder[count].transform(input_data[i]))
count=count+1
input_data_encoded=np.array(input_data_encoded)
print('Predicted traffic: ', int(regressor.predict(input_data_encoded)[0])) # 打印分类结果
五、无监督学习--聚类：
无监督学习是一种对不含标记的数据建立模型的机器学习方法。最常见的无监督学习方法是聚类，当需要把无标记的数据分成几种集群时，就需要用来分析，这些集群通常是根据某种相似度指标进行划分，如欧式距离。无监督学习广泛应用于各种领域，如数据挖掘、医学影像、股票市场分析、计算机视觉、市场细分等。
1. k-means算法聚类数据：
这个算法常常利用数据的不同属性将输入数据划分为k组，分组使用最优化的技术实现，即让组内数据点与改组中心点的距离平方和最小。
import numpy as np
import matplotlib.pyplot as plt
from sklearn import metrics
from sklearn.cluster import KMeans
import utilities
加载数据，定义集群的数量。
data=utilities.load_data('data_multivar.txt')
num_clusters=4
显示数据：
plt.figure()
plt.scatter(data[:,0], data[:,1], marker='o', facecolors='none', edgecolors='k', s=30)
x_min, x_max=min(data[:,0])-1, max(data[:,0])+1
y_min, y_max=min(data[:,1])-1, max(data[:,1])+1
plt.title('Input data')
plt.xlim(x_min, x_max)
plt.ylim(y-min, y_max)
plt.xticks(())
plt.yticks(())
训练模型：
kmeans=KMeans(init='k-means++', n_clusters=num_clusters, n_init=10)
kmeans.fit(data)
数据训练之后需要可视化边界：
step_size=0.01
x_min, x_max=min(data[:,0])-1, max(data[:,0])+1
y_min, y_max=min(data[:,1])-1, max(data[:,1])+1
x_values, y_values=np.meshgrid(np.arange(x_min, x_max, step_size), np.arange(y_min, y_max, step_size))
predicted_labels=kmeans.predict(np.c_[x_values.ravel(), y_values.ravel()])
绘出结果：
predicted_labels=predicted_labels.reshape(x_values.shape)
plt.figure()
plt.clf()
plt.imshow(predicted_labels, interpolation='nearest', extent=(x_values.min(), x_values.max(), y_values.min(), y_values.max()), cmap=plt.cm.Paired, aspect='auto', origin='lower')
plt.scatter(data[:,0], data[:,1], marker='o', facecolors='none', edgecolors='k', s=30)
把中心点画在图形上：
centroids=kmeans.cluster_centers_
plt.scatter(centroids[:,0], centroids[:,1], marker='o', s=200, linewidths=3, color='k', zorder=10, facecolors='black')
x_min, x_max=min(data[:,0])-1, max(data[:,0])+1
y_min, y_max=min(data[:,1])-1, max(data[:,1])+1
plt.title('Centoids and boundaries obtained using KMeans')
plt.xlim(x_min, x_max)
plt.ylim(y_min, y_max)
plt.xticks(())
plt.yticks(())
plt.show()
2. 用矢量量化压缩图片：
k-means聚类的主要应用之一就是矢量量化，也就是舍入的N维版本。矢量量化用于图片压缩，就可以使用比原始图像少的比特数来存储像素。
import argparse
import numpy as np
from scipy import misc
from sklearn import cluster
import matplotlib.pyplot as plt
一个自定义的函数，用于解析输入数据，把图片和每个像素被压缩的比特数传入：
def build_arg_parser():
parse=argparse.ArgumentParse(description='Compress the input image using clustering')
parse.add_argument('--input-file', dest='input_file', required=True, help='Input image')
parse.add_argument('--num-bits', dest='num_bits', required=False, type=int, helf='Number of bits used to repersent each pixel')
return parse
创建自定义函数，用来压缩输入图片：
def compress_image(img, num_clusters):
X=img.reshape((-1,1)) # 将输入图片转换成（样本量，特征量）数组
kmeans=cluster.KMeans(n_clusters=num_clusters, n_init=4, random_state=5)
kmeans.fit(X)
centroids=kmeans.cluster_centers_.squeeze()
labels=kmeans.labels_
# 为每个数据配置离它最近的中心点，并转变为图片的形状
input_image_compressed=np.choose(labels, centroids).reshape(img.shape)
return input_image_compressed
压缩完图片后，显示图片以观察图片质量：
def plot_image(img, title):
vmin=img.min()
vmax=img.max()
plt.figure()
plt.title(title)
plt.imshow(img, cmap=plt.cm.gray, vmin, vmax=vmax)
定义主函数，将输入数据传入并进行处理，然后提取输出图片：
if __name__=='__main__':
args=build_arg_parse().parse_args()
input_fiel=args.input_file
num_bits=args.num_bits
if not 1<=num_bits<=8:
raise TypeError('Number of bits should be between 1 and 8')
num_clusters=np.power(2, num_bits)
compression_rate=round(100*(8.0-args.num_bits)/8.0, 2) # 计算压缩率
print('\nThe size of the image will be reduced by a factor of ',8.0/args.num_bits)
print('\nCompression rate='+str(compression_rate)+'%')
加载输入图片：
input_image=mosc.imread(input_file, True).astype(np.uint8)
plot_image(input_image, 'Original image')
用输入参数压缩图片：
input_image_compressed=compress_image(input_image, num_clusters)
plot_image (input_image_compressed, 'Compressed image; compression rate =' +str(compression_rate) +'5')
plt.show()
3. 均值漂移聚类模型：
该算法把数据点的分布看成是概率密度函数，希望在特征空间中根据函数分布特征找出数据点的模式。这些模式就对应一群群局部最密集分布的点。优点是不需要事先确定集群的数量。
首先加载数据：
X=utilities.load_data('data_nultivar.txt')
通过指定参数创建一个均值漂移模型：
bandwidth=estimate_bandwidth(X, quantile=0.1, n_samples=len(X))
meanshift_estimator=MeanShift(bandwidth=bandwidth, bin_seeding=True)
训练模型：
meanshift_estimator.fit(X)
提取标记：
labels=meanshift_estimator.labels_
从模型中提取集群的中心点，然后打印集群数量：
centroids=meanshift_estimator.cluster_centers_
num_clusters=len(np.unique(labels))
print('Number of clusters in input data= ', num_clusters)
把集群可视化：
import matplotlib.pyplot as plt
from intertools import cycle
plt.figure()
markers='.*xv'
迭代数据点并画出：
for i,marker in zip(range(num_clusters), markers):
plt.scatter(X([labels==i, 0], X[labels==i, 1], marker=marker, color='k')
centriod=centroids[i]
plt.plot(centroid[0], centroid[1], marker='o', markerfacecolor='k', markeredgecolor='k', markersize=15)
plt.title('Clusters and their centroids')
plt.show()
4. 用凝聚层次聚类进行数据分组：
层次聚类hierarchical clustering是一组聚类算法，通过不断地分解或合并集群来构建树状集群。层次聚类可以用一棵树表示。
层次聚类可以是自下而上的，也可以是自上而下的。在自下而上的算法中，每个数据点都被看作是一个单独的集群，这些集群不断合并，直到所有的集群都合并成一个巨型集群，这被称为凝聚层次聚类。与之相反，自上而下的层次算法是从一个巨大集群开始，不断分解，直到所有的集群变成一个单独的数据点。
import numpy as np
import matplotlib.pyplot as plt
from sklearn.cluster import AgglomerativeClustering
from sklearn.neighbors import kneighbors_graph
定义一个实现凝聚层次聚类的函数：
def perform_clustering(X, connectivity, title, num_clusters=3, linkage='ward'):
plt.figure()
model = AgglomerativeClustering (linkage=linkage, connectivity=connectivity, m_clusters=num_clusters)
提取标记，然后指定不同聚类在图像中的标记：
labels=model.labels_
markers='.vx'
迭代数据，用不同的标记把聚类的点画在图像中：
for i,marker in zip(range(num_clusters), markers):
plt.scatter(X[labels==i,0],X[labels==i,1], s=50, marker=marker, color='k', facecolors='none')
plt.title(title)
为了演示凝聚层次聚类的优势，对一些在空间中是否连接在一起、但彼此却非常接近的数据进行聚类。希望连接在一起的数据可以聚成一类，而不是在空间上非常接近的点聚成一类。定义一个函数来获取一组呈螺旋状的数据点：
def get_spiral(t, noise_amplitude=0.5):
r=t;x=r*np.cos(t);y=r*np.sin(t)
return add_noise(x,y,noise_amplitude)
在上述函数中增加了一些噪声，可以增加一些不确定性。噪声函数为：
def add_noise(x,y,amplitude):
X=np.concatenate((x))
X+=amplitude*np.random(2, X.shape[1])
return X.T
再定义一个函数来获取位于玫瑰曲线上的数据点：
def get_rose(t, noise_amplitude=0.02):
k=5;r=np.cos(k*t)+0.25;x=r*np.cos(t);y=r*np.sin(t)
return add_noise(x,y,noise_amplitude)
为了增加多样性，再定义一个hypotrochoid函数：
def get_hypotrochoid(t, noise_amplitude=0):
a,b,h=10.0,2.0,4.0
x=(a-b)*np.cos(t)+h*cos((a-b)/b*t)
y=(a-b)*np.sin(t)-h*sin((a-b)/b*t)
return add_noise(x,y,0)
定义主函数：
if __name__=='__main__':
n_samples=500 # 生成样本数据
np.random.seed(2)
t=2.5*np.pi*(1+2*np.random.rand(1, n_samples))
X=get_spiral(t)
connectivity=None # 不考虑螺旋形的数据连续性
perform_clustering(X, connectivity, 'No connectivity')
connectivity=kneighbors_graph(X, 10, include_self=False)
perform_clustering(X, connectivity, 'K-Neighbors connectivity')
plt.show()
使用连接特性，可以让连接在一起的数据合成一组，而不是按照它们在螺旋线上的位置进行聚类。
5. 评价聚类算法的聚类效果：
度量聚类算法的效果采用轮廓系数Silhouette Coefficient，定义为：(x-y)/max(x,y)。其中，x表示在同一个集群中某个数据点与其他数据点的平均距离，y表示某个数据点与最近的另一个集群的所有点的平均距离。
加载数据：
data=utilities.load_data('data_perf.txt')
为了确定集群的最佳数量，迭代一系列的值，找出其中的峰值：
scores=[]
range_values=np.arange(2,10)
for i in range_values:
kmeans=KMeans(init='k-means++', n_clusters=i, n_init=10)
kmeans.fit(data)
score=metrics.silhouette_score(data, kmeans.labels_, metric='euclidean', sample_size=len(data))
print('\nNumber of clusters=', i)
print('Silhouette score', score)
scores.append(score)
画出图形并找出峰值：
plt.figure()
plt.bar(range_values, scores, width=0.6, color='k', align='center') #画得分条形图
plt.title('Silhouette score vs number of clusters')
plt.figure()
plt.scatter(data[:,0], data[:,1], color='k', s=30, marker='o', facecolors='none')
x_min,x_max=min(data[:,0])-1,max(data[:,0])+1
y_min,y_max=min(data[:,1])-1,max(data[:,1])+1
plt.title('Input data')
plt.xlim(x_min, x_max)
plt.ylim(y_min, y_max)
plt.xticks(())
plt.yticks(())
plt.show()
6. 用DBSCAN算法指定估算集群数量：
带噪声的基于密度的聚类方法DBSCAN（Density-Based Spatial Clustering of Applications with Noise）将数据点看成是紧密集群的若干组。如果某个点属于一个集群，那么就应该有许多点也属于同一个集群。该方法中有一个spsilon参数，可以控制这个点到其他点的最大距离，如果两个点的距离超过了参数epsilon的值，它们就不可能在一个集群。此方法的主要优点是可以处理异常点，如果一些点位于数据稀疏区域，就会把这些点最为异常点，而不会强制放入一个集群中。
加载数据：
input_file='data_perf.txt'
X=load_data(input_file)
需要找到最佳集群数量参数，先初始化一些变量：
eps_grid=np.linspace(0.3, 1.2, num=10)
silhouette_scores=[]
eps_best=eps_grid[0]
silhouette_score_max=-1
model_best=None
labels_best=None
搜索参数空间：
for eps in eps_grid:
model=DBSCAN(eps=eps, min_samples=5).fot(X) # 训练DBSCAN聚类模型
labels=model.labels_ # 提取标记
每次迭代都需要提取性能指标：
silhouette_score=round(metrics.silhouette_score(X, labels), 4)
silhouette_scores.append(silhouette_score)
print('Epsilon: ', eps, '--> silhouette score: ',silhouette_score)
需要保存指标的最佳得分和对应的epsilon值：
if silhouette_score>silhouette_score_max:
silhouette_score_max=silhouette_score
eps_best=eps
model_best=model
labels_best=labels
画出条形图：
plt.figure()
plt.bar(eps_grid, silhouette_scores, width=0.05, color='k', align='center')
plt.title('Silhouette score vs spsilon')
print('\nBest epsilon=', eps_best) # 打印最优参数
把最优的模型和标记保存：
model=model_best
labels=labels_best
有些数据点还没有分配集群，需要识别：
offset=0
if -1 in labels:
offset=1
提取集群的数量：
num_clusters=len(set(lables))-offset
print(]\nEstimated number of clusters=', num_clusters)
提取核心样本：
mask_core=np.ziros(labels.shape, dtype=np.bool)
mask_core[model.core_sample_indices_]=True
接下来将集群结果可视化，首先提取独特的标记集合，然后分配不同的标记：
plt.figure()
labels_uniq=set(labels)
markers=cycle('vo^s<>')
用迭代法把每个集群的数据点用不同的标记画出来：
for cur_label,marker in zip(labels_uniq, markers):
if cur_label==-1:
marker='.' # 用黑点表示未分配的数据点
cur_mask=(labels==cur_label)
cur_data=X[cur_mask & mask_core]
plt.scatter(cur_data[:,0], cur_data[:,1], marker=marker, edgecolors='black', s=96, facecolors='none')
cur_data=X[cur_mask & ~mask_core]
plt.scatter(cur_data[:,0], cur_data[:,1], marker=marker, edgecolors='black', s=32)
plt.title('Data separated into clusters')
plt.show()
7. 近邻传播聚类算法：
近邻传播聚类Affinity Propagation算法会找出数据中每个集群的代表性数据点，会找到数据点间的相似度度量值，并把所有数据点看出潜在的代表性数据点，也称为取样器exemplar。
用这种方法分析骨架的波动找出公司行为的相似性。
import json
import datetime
import numpy as np
import matplotlib.pyplot as plt
from sklearn import covariance,cluster
from matplotlib.finance import quotes_historical_yahoo_ochl as quotes_yahoo
加载数据：
symbol_file='symbol_map.json'
从符号映射文件中读取数据：
with open(symbol_file, 'r') as f:
symbol_dict=json.loads(f.read())
symbols,names=np.array(list(symbol_dict.items())).T
指定分析的时间段，作为输入数据的起止时间：
start_date=datetime.datetime(2004,4,5)
end_date=datetime.datetime(2007,6,2)
读取输入数据：
quotes=[quotes_yahoo(symbol, start_date, end_date, asobject=True) for symbol in symbols]
使用每天的开盘价和收盘价的差异来分析数据：
opening_quotes=np.array([quote.open for quote in quotes]).astype(np.float)
closing_quotes=np.array([quote.close for quote in quotes]).astype(np.float)
delta_model=covariance.GraphLassoCV() # 计算每日股价波动
建立一些协方差图模型：
edge_model=covariance.GraphLassoCV()
在使用数据前对它进行标准化：
X=delta_quotes.copy().T
X/=X.std(axis=0)
用数据训练模型：
with np.errstate(invalid='ignore'):
edge_model.fit(X)
建立聚类模型：
_, labels=clusters.affinity_propagation(edge_model.covariance_)
num_labels=labels.max()
打印结果：
for i in range(num_labels+1):
print('Cluster', i+1, '-->', ', '.join(names[labels==i])
8. 客户细分模型：
无监督学习的主要应用之一是市场细分，这对广告投放、库存管理、配送策略等市场行为都非常有用。
import csv
import numpy as np
import matplotlib.pyplot as plt
from sklearn.cluster import MeanShift, estimate_bandwidth
载入数据：
input_file='wholesale.csv'
file_reader=csv.reader(open(input_file, 'rb'), delimiter=',')
X=[]
for count,row in enumerate(file_reader):
if not count:
names=row[2:]
continue
X.append([float(x) for x in row[2:]])
X=np.array(X)
建立一个均值漂移聚类模型：
bandwidth=estimate_bandwidth(X, quantile=0.8, n_samples=len(X))
meanshift_estimator=MeanShift(bandwidth=bandwidth, bin_seeding=True)
meanshift_estinator.fit(X)
labels=meanshift_estimator.labels_
centroids=meanshift_estimator.cluster_centers_
num_clusters=len(np.unique(labels))
print('\nNumber of clusters in input data=', num_clusters)
打印获得的集群中心：
print('\nCentroids of clusters: ')
print('\t'.join([name[:3] for name in names])
for centroid in centroids:
print('\t'.join([str(int(x)) for x in centroid])
把两个特征milk和groceries的聚类结果可视化，以获取直观的输出：
centroids_milk_groceries=centroids[:,1:3]
plt.figure()
plt.scatter(centroids_milk_groceries[:,0], centroids_nilk_groceries[:,1], s=100, edgecolors='k', facecolors='none')
offset=0.2
plt.xlim(centroids_milk_groceries[:,0].min()-offset*centroids_milk_groceries[:,0].ptp(), centroids_milk_groceries[:,0].max()-offset*centroids_milk_groceries[:,0].ptp(),)
plt.ylim(centroids_milk_groceries[:,1].min()-offset*centroids_milk_groceries[:,1].ptp(), centroids_milk_groceries[:,1].max()-offset*centroids_milk_groceries[:,1].ptp(),)
plt.title('Centroids of clusters for nilk and groceries')
plt.show()
六、构建推荐引擎：
推荐引擎是一个预测用户兴趣点的模型。通过预测当前用户可能会喜欢的内容，将相应的物品从数据库中筛选出来，这样的推荐引擎可以有助于将用户和数据集中的合适内容连接起来。如果有一个庞大的商品目录，而用户不可能查找所有的相关内容，通过推荐合适的内容，可以增加用户消费。
推荐引擎通常用协同过滤collaborative filtering或基于内容的过滤content-based filtering来产生一组推荐。协同过滤从当前用户过去的行为和其他用户对当前用户的评分来构建模型，然后使用这个模型来预测这个用户可能感兴趣的内容；而基于内容的过滤用商品本身的特征来给用户推荐更多商品，商品间的相似度是模型主要的关注点。
1. 构建机器学习流水线：
scikit-learn库中包含了构建机器学习流水线的方法，只需要指定函数来构建一个组合对象，使数据通过整个流水线。这个流水线包括预处理、特征选择、监督式学习、非监督式学习等函数。
from sklearn.datasets import samples_generator
from sklearn.ensemble import RandomForestClassifier
from sklearn.feature_selection import SelectKBest,f_regression
from sklrean.pipeline import Pipeline
生成示例数据：
X,y=samples_generator.make_classification(n_informative=4, n_features=20, n_redundant=0, random_state=5)
代码生成了一个20维的特征向量，这是一个默认值。可以通过修改n_features参数来修改特征向量的维数。
在数据点进一步操作之前选择k个最好的特征，设置k的值为10：
selector_k_best=SelectKBest(f_regression, k=10)
用随机森林分类器分类数据：
classifier=RandomForestClassifier(n_estimators=50, max_depth=4)
创建流水线：
pipeline_classifier=Pipeline(['selector', selector_k_best), ('rf', classifier)])
可以在流水线中为模块指定名称。还可以更新参数：
pipeline_classifier.set_params(selector__k=6, rf__n_estimators=25)
训练分类器：
pipeline_classifier.fit(X,y)
预测输出结果：
prediction=pipeline_classifier.predict(X)
print('\nTredictions:\n', prediction)
评价分类器的性能：
print('\nScore: ', pipeline_classifier.score(X,y))
还可以查看哪些特征被选中，并输出显示：
features_status=pipeline_classifier.named_steps['selector'].get_support()
selected_features=[]
for count,item in enumerate(features_status):
if item:
selected_features.append(count)
print('\nSelected features (0-indexed):', ', ', '.join([str(x) for x in selected_features])
选择k个最好的特征，可以处理较小维度的数据，这对减小计算复杂度很有用。选择k个最佳特征的方式是基于单变量的特征选择，选择过程是先进行单变量统计测试，然后从特征向量中抽取最优秀的特征。做了这些测试后，向量空间中的每个特征酱油一个评价分数。基于这些评价分数，选择最好的k个特征。一旦抽取出k个特征，一个k位的特征向量就形成了，可以将这个特征向量用于随机森林分类器的输入训练数据。
2. 寻找最近邻：
最近邻模型是指一个通用算法类，目的是根据训练数据集中的最近邻数量来做决策。
输入数据：
X=np.array([[1,1],[1,3],[2,2],[2.5,5],[3,1],[4,2],[2,3.5],[3,3],[3.5,4]])
目标是对于任意给定点找到3个最近邻：
num_neighbors=3
定义一个随机数据点，这个数据点不包括在输入数据中：
input_point=[2.6,1.7]
画出数据分布：
plt.figure()
plt.scatter(X[:,0], X[:,1], marker='o', s=25, color='k')
为了寻找最近邻，用适合的参数定义一个NearestNeighbors对象，并用输入数据训练该对象：
knn=NearstNeighbors(n_neighbors=num_neighbors, algorithm='ball_tree').fit(X)
计算输入点与输入数据中所有点的距离：
distances,indices=knn.kneighbors(input_point)
打印出k个最邻近：
for rank,index in enumerate(indices[0][:num_neighbors]):
print(str(rank+1)+'-->', X[index])
indices数组是一个已排序数组，仅需要解析并打印出数据点。
画出输入数据点，并突出显示k个最近邻：
plt.figure()
plt.scatter(X[:,0], X[:,1], marker='o', s=25, color='k')
plt.scatter(X[indices][0][:][:,0], X[indices][0][:][:,1], marker='o', s=150, color='k', facecolors='none')
plt.scatter(input_point[0], input_point[1], marker='x', s=150, color='k', facecolors='none')
plt.show()
3. 构建一个KNN分类器：
KNN（k-nearest neighbors）是用k个最近邻的训练数据集来寻找未知对象分类的一种算法。如果希望找到未知数据点属于哪个类，可以找到KNN并做一个多数表决。
import numpy as np
import matplotlib.pyplot as plt
import matplorlib.cm as cm
form sklearn import neighbors,dataset
from utilities import load_data
以data.txt文件为输入数据，加载：
input_file='data.txt'
data.load_data(input_file)
X,y=data[:,:-1],data[:,-1].astype(np.int)
前两列包含输入数据，最后一列包含标签，因此分别将其用X和y两个变量表示。
将输入数据可视化：
plt.figure()
plt.title('Input datapoint')
markers='^sov<>hp'
mapper=np.array([markers[i] for i in y])
for i in range(X.shape[0]):
plt.scatter(X[i,0],X[i,1],marker=mapper[i], s=50, edgecolors='black', facecolors='none')
迭代所有的数据点并用适合的标记区分不同的类。
为了构建分类器，需要指定考虑的最近邻的个数：
num_neighbors=10
为了将边界可视化，需要定义一个网格，用这个网格评价该分类器：
h=0.01 # 定义网格步长
创建KNN分类器并进行训练：
classifier=neighbors.KNeighborsClassifier(num_neighbors, weights='distance')
classfifier.fit(X,y)
生成一个网格来画出边界：
x_min,x_max=X[:,0].min()-1, X[:,0].max()+1
y_min,y_max=X[:,1].min()-1, X[:,1].max()+1
x_grid,y_grid=np.meshgrid(np.arange(x_min, x_max, h)
np.arange(y_min, y_max, h)
评价分类器对所有点的输出：
predicted_values=classifier.predict(np.c_[x_grid.ravel(), y_grid.ravel()])
将其画出：
predicted_values=predicted_values.reshape(x_grid.shape)
plt.figure()
plt.pcolormesh(x_grid, y_grid, predicted_values, cmap=cm.Pastell)
将训练数据点覆盖在其上：
for i in range(X.shape[0]):
plt.scatter(X[i,0], X[i,1], marker=mapper[i], s=50, edgecolors='black', facecolors='none')
plt.xlim(x_grid.min(), x_grid.max())
plt.ylim(y_grid.min(), y_grid.max())
plt.title('k nearest neighbors classifier boundaries')
测试数据点，查看分类器能否准确分类，并画出：
test_datapoint=[4.5,3.6]
plt.figure()
plt.title('Test datapoint')
for i in range(X.shape[0]):
plt.scatter(X[i,0], X[i,1], marker=lapper[i], s=50, edgecolors='black', facecolors='none')
plt.scatter(test_datapoint[0], test_datapoint[1], marker='x', linewidth=3, s=200, facecolors='black')
用以下模型提取KNN：
dist, indices=classifier.kneighbors(test_datapoint)
画出KNN并突出显示：
plt.fugure()
plt.title('k nearest neighbors')
for i in indices:
plt.scatter(X[i,0], X[i,1], marker='o', linewidth=3, s=100, facecolors='black')
plt.scatter(test_datapoint[0], test_datapoint[1], marker='x', linewidth=3, s=200, facecolors='black')
for i in range(X.shape[0]):
plt.scatter(X[i,0], X[i,1], marker=mapper[i], s=50, edgecolors='black', facecolors='none')
plt.show()
输出分类器输出结果：
print('Predicted output: ', classifier.predict(test_datapoint)[0]
KNN分类器存储了所有可用的数据点，并根据相似度指标来对新的数据点进行分类，这个相似度指标通常以距离函数的形式度量。该算法是一个非参数化技术，也就是说它在进行计算前并不需要找出任何隐含的参数，只需要选择k的值。
一旦找出KNN，就会做一个多数表决。一个新数据点通过KNN的多数表决来进行分类，这个数据点会被分到KNN最常见的类中。如果将k的值设置为1，就变成一个最近邻分类器。
4. 构建一个KNN回归器：
import numpy as np
import matplotlib.pyplot as plt
form sklearn import neighbors
生成一些服从正态分布的样本数据：
amplitude=10
num_points=100
X=amplitude*np.random.rand(num_points,1)-0.5*amplitude
在数据点中加入一些噪声引入一些随机性：
y=np.sinc(X).ravel()
y+=0.2*(0.5-np.random.rand(y.size))
数据可视化：
plt.figure()
plt.scatter(X, y, s=40, c='k', facecolors='none')
plt.title('Input data')
定义更密集的网格：
x_values=np.linespace(-0.5*amplitude, 0.5*amplitude, 10*num_points)[:, np.newaxis]
定义最近邻的个数：
n_neighbors=8
用前面定义的参数初始化并训练KNN回归器：
knn_regressor=neighbors.KNeighborsRegressor(n_neighbors, weights='diatance')
y_values=knn_regressor.fit(X, y).predict(x_values)
将输入数据和输出数据交叠在一起，以查看回归器的性能：
plt.figure()
plt.scatter(X, y, s=40, c='k', facecolors='none', label='input data')
plt.plot(x_values, y_values, c='k', linestyle='--', label='predicted values')
plt.xlim(X.min()-1, X.max()+1)
plt.ylim(y.min()-0.2, y.max()+0.2)
plt.aopxis('tight')
plt.legend()
plt.title('K Nearest Neighbors Regressor')
plt.show()
回归器的目标是预测连续值的输出。示例中用到sinc函数来演示KNN回归器，该函数也称基本正弦函数，定义为sinc(x)=sin(x)/x。当x=0时，变成0/0不定式，因此需要计算函数在x无限趋近于0时函数的极限。
5. 计算欧式距离分数：
为了构建一个推荐引擎，需要定义相似度指标，以便找到与数据库中特定用户相似的用户。可以计算两个数据点之间的欧几里得距离。定义用于计算欧几里得分数的函数：
def euclidean_scire(dataset, user1, user2):
if user1 not in dataset:
raise TypeError('User '+user1+' not present in the dataset')
if user2 not in dataset:
raise TypeError('User '+user2+' not present in the dataset')
# 提取两个用户均评过分的电影
rated_by_both={}
for item in dataset[user1]:
if item in dataset[user2]:
rated_by_both[item]=1
# 如果没有两个用户共同评过分的电影，说明两个用户之间没有相似度
if len(rated_by_both)==0:
return 0
# 对于每个共同评分，只计算平方和和平方根，并将该值归一化到0与1之间
squared_differences=[]
for item in dataset[user1]:
if item in dataset[user2]:
squared_differences.append (np.square(dataset[user1][item] - dataset[user2][item]))
return 1/(1+np.sqrt(np.sum(aquared_differences)))
如果评分相似，那么平方和的差别就会很小，因此评分就会很高。
加载数据文件：
if __name__=='__main__':
data_file='movie.json'
with open(data_file,'r') as f:
data=json.loads(f.read())
# 假定两个随机用户，计算其欧式距离分数
user1='John Carson'
user2='Michelle Peterson'
print(euclidean_score(data, user1, user2)
6. 计算皮尔逊相关系数：
定义计算两个用户之间的皮尔逊相关度系数的函数：
def pearson_scire(dataset, user1, user2):
if user1 not in dataset:
raise TypeError('User '+user1+' not present in the dataset')
if user2 not in dataset:
raise TypeError('User '+user2+' not present in the dataset')
# 提取两个用户均评过分的电影
rated_by_both={}
for item in dataset[user1]:
if item in dataset[user2]:
rated_by_both[item]=1
num_ratings=len(rated_by_both)
# 如果没有两个用户共同评过分的电影，说明两个用户之间没有相似度
if num_ratings==0:
return 0
# 对于每个共同评分，只计算平方和和平方根，并将该值归一化到0与1之间
user1_squared_sum=np.sum([np.square(dataset[user1][item]) for item in rated_by_both])
user2_squared_sum=np.sum([np.square(dataset[user2][item]) for item in rated_by_both])
# 计算数据集的乘积之和
product_sum=np.sum([dataset[user1][item]*dataset[user2][item] for item in rated_by_both])
# 计算皮尔逊相关系数需要的各种元素：
Sxy=product_sum-(user1_sum*user2_sum/num_ratings)
Sxx=user1_squared_sum-np.square(user1_sum)/num_ratings
Syy=user2_squared_sum-np.square(user2_sum)/num_ratings
if Sxx*Syy==0:
return 0
return Sxy/np.sqrt(Sxx*Syy)
定义main函数并计算两个用户之间的皮尔逊相关系数：
if __name__=='__main__':
data_file='movie.json'
with open(data_file,'r') as f:
data=json.loads(f.read())
# 假定两个随机用户，计算其欧式距离分数
user1='John Carson'
user2='Michelle Peterson'
print(pearson_score(data, user1, user2)
7. 寻找数据集中的相似用户：
构建推荐引擎中一个非常重要的任务是寻找相似的用户，即为某位用户生成推荐信息可以同时推荐给与其相似的用户。
from pearson_score import pearson_score
定义一个函数，用于寻找与输入用户相似的用户，该函数有3个输入参数，数据库、输入用户、寻找相似用户的个数。
def find_similar_users(dataset, user, num_users):
if user not in dataset:
raise TypeError('User '+user+' not present in the dataset')
# 计算所有用户的皮尔逊相关度
scores=np.arrsy([x, pearson_score(dataset, user, x)] for in dataset if user != x])
scores_sorted=np.argsort(scores[:,1]) # 评分按照第二列排序
scored_sorted_dec=scores_sorted[::-1] # 评分按照降序排列
# 提出k个最高分
top_k=scored_sorted_dec[0:num_users]
return scores[top_k]
定义main函数，加载输入数据库：
if __name__='__main__':
data_file='movie.json'
with open(data_file, 'r') as f:
data=json.loads(f.read())
user='John Carson'
print('\nUsers similar to '+user+':\n'
similar_users=find_similar_users(data, user, 3) # 查找3个相似用户
print('User\t\t\tSimilarity scores\n'
for item similar_users:
print(item[0], '\t\t', round(float(item[1]), 2)
8. 生成电影推荐：
import json
import numpy as np
from euclidean_score import euclidean_score
from pearson_score omport pearson_score
from find_similar_users import find_similar_users
定义一个给用户生成推荐的函数：
def generate_recommendations(dataset, user)：
if user nor in dataset:
raise TypeError('User '+user+' not present in the dataset')
total_scores={}
similarity_sums={}
for u in [x for x in dataset if x!=user]:
similarity_score=pearson_score(dataset, user, u)
if similarity_score<=0:
continue
# 找到还未被该用户评分的电影
for item in [x for x in dataset[u] if x not in dataset[user] or dataset[user][x]==0]:
total_scores.update({item:dataset[u][item]*similarity_score})
similarity_sums.update({item:similarity_score})
# 如果该用户看过数据库中的所有电影，就不为用户推荐
if len(total_scores)==0:
return ['No recommendations possible']
# 生成一个电影评分标准化列表
movie_ranks=np.array([total/similarity_sums[item], item]) for item,total in total_scores.items()])
# 根据皮尔逊相关系数进行降序排列：
movie_ranks=movie_ranks[np.argsort(movie_ranks[:,0])[::-1]]
# 提取推荐电影
recommendations=[movie for _,movie in movie_ranks]
return recommendations
定义主函数，加载数据集：
if __name__=='__main__':
data_file='movie_ratings.json'
with open(data_file,'r') as f:
data=json.loads(f.read())
user='Michael Henry'
print('\nRecommendations for '+user+':')
movies=generate_recommendations(data,user)
for i,movie in enumerate(movies):
print(str(i+1)+'. '+movie)
user='John Carson'
print('\nRecommendations for '+user+':')
movies=generate_recommendations(data,user)
for i,movie in enumerate(movies):
print(str(i+1)+'. '+movie)
七、分析文本数据：
自然语言处理NLP（Natural Language Processing）是现代人工智能的重要领域，最用于搜索引擎、情感分析、主题建模、词性标注、实体识别等。
Python中的NLTK（Natural Language Toolkit）包可用于基于机器学习的NLP检测文本数据。
1. 用标记解析的方法预处理数据：
标记解析是将文本分割成一组有意义的片段的过程，这些片段被称作标记。可以将一段文字分割成单词或者句子，可以自定义将输入文本分割成有意义的标记。
from nltk.tokenize import sent_tokenize
sent_tokenize_list=sent_tokenize(text)
NLTK附带了几个不同的单词解析器：
from nltk.tokenize import word_tokenize
word_list=word_tokenize(text)
NLTK可以使用PunktWord单词解析器，以标点符号分割文本：
from nltk.tokenize import PunktWordTokenizer
punkt_tokenizer=PinktWorldTokenizer()
print(punkt_tokenizer.tokenize(text))
如果需要将标点符号保留到不同的句子中，用WordPunct标记解析器：
from nlyk.tokenize import WordPunctTokenizer
wordpunct_tokenizer=WordPunctTokenizer()
print(wordpunct_tokenizer.tokenize(text))
2. 提取文本数据的词干：
处理文本文档时，可能会碰到单词的不同形式。在文本分析中，提取这些单词的原形很有用。
3. 用分块的方法划分文本：
分块是指基于任意随机条件将输入文本分割成块，当处理非常大的文本文档时就需要分块，以便下一步分析。
4. 创建词袋模型：
如果需要处理包含数百万单词的文本文档，需要将其转换成某种数值表示形式，以便让机器用这些数据来学习算法。这些算法需要数值数据，以便可以对这些数据进行分析，并输出有用的信息。这里需要用到词袋（bag-of-words）。词袋是从所有文档的所有单词中学习词汇的模型。学习之后，词袋通过构建文档中所有单词的直方图来对每篇文档进行建模。
5. 创建文本分类器：
文本分类的目的是将文本文档分为不同的类，这是NLP中非常重要的分析手段。这里将使用tf-idf统计数据，表示词频-逆文档频率。
tf-idf技术常用于信息检索领域，目的是了解文档中每个单词的重要性。词频表示一个单词在给定文档中出现的频次，为了规范化，用频次除以文档中所有单词的个数。逆文档频率表示给定单词的重要性，为了抗衡那些经常出现的单词频率，需要用一个系数将其权重变小，要计算文档的总数目除以该单词出现的文档数目的比值，逆文档频率对该比值取对数值。
6. 识别性别：
NLP中，通过姓名识别性别是一个有趣的任务。
7. 分析句子的情感：
情感分析是NLP最受欢迎的应用之一，是指确定一段给定的文本是积极还是消极，还可以使用中性选项。
8. 用主题建模识别文本的模式：
主题建模是指识别文本数据隐藏模式的过程，目的是发现一组文档的隐藏主题结构。
八、语音识别：
语音识别是指识别和理解口语的过程。输入音频数据，语音识别器将处理这些数据，从中提取出有用的信息。语音识别有很多实际应用，如声音控制设备等。
1. 读取和绘制音频数据：
实际的音频信号是复杂的连续波形，为了将其保存为数字化形式，需要对音频信号进行采样并将其转换成数字。
import numpy as np
import matplotlib.pyplot as plt
from scopy.io import wavfile
使用wavfile包读取音频文件：
samping_freq,audio=wavfile.read('input_read.wav')
输出相关参数：
print('\nShape: ', audio.shape)
print('Datatype: ', audo.dtype)
print('Duration: ', round(audio.shape[0]/float(sample_freq), 3), 'seconds')
该音频信号存储为16位有符号整型数据，标准化这些值：
audio=audio/(2.**15)
提取前30个值，并绘图：
audio=audio[:30]
X轴为时间轴，且按采样频率因子进行缩放：
x_values=np.arrange(0,len(audio),1)/float(sampling_freq)
将单位转换为秒：
x_values*=1000
绘制图形：
plt.plot(x_values, audio, color='black')
plt.xlabel('Time (ms)')
plt.ylabel('Amplitude')
plt.title('Audio signal')
plt.show()
2. 将音频信号转换为频域：
音频信号可以认为是不同频率、幅度和相位的正弦波的复杂混合。可以使用傅里叶变换转换到频域。
import numpy as np
from scipy.io import wavfile
import matplotlib.pyplot as plt
读取音频文件：
sampling_freq,audio=wavfile.read('input_freq.wav')
对信号进行标准化：
audio=audio/(2.**15)
音频信号就是一个NumPy数组，提取其长度：
len_audio=len(audio)
傅里叶变换是关于中心点对称的，因此只需要转换信号的前半部分。因为要提取信号的概率，需要要将信号的值平方：
transformed_signal=np.fft.fft(audio)
half_length=np.ceil((len_audio+1)/2.0)
transformed_signal=abs(transformed_signal[0:half_length])
transformed_signal/=float(len_audio)
transformed_signal**=2
提取信号的长度：
len_ts=len(transformed_signal)
根据信号的长度将信号乘以2：
if len_audio%2:
transformed_signal[1:len_ts]*=2
else:
transformed_signal[1:len_ts-1]*=2
功率信号用下式求得：
power=10*np.log10(transformed_signal)
X轴是时间轴。接下来需要根据采样频率对其进行缩放，并将其转换成秒：
x_values=np.arange(0,half_length,1)*(sampling_freq/len_audio)/1000.0
绘制该信号：
plt.figure()
plt.plot(x_values, power, color='black')
plt.xlabel('Freq(in kHz)')
plt.ylabel('Power(in dB)')
plt.show()
3. 自定义参数生成音频信号：
import numpy as np
from scipy.io import wavfile
import matplotlib.pyplot as plt
定义输出文件，用于存储生成的音频：
output_file='output_gen.wav'
指定音频生成参数，希望生成一个长度3s的信号，采样频率44100Hz，音频频率587Hz，时间轴的值从-2xpi到2xpi：
duration=3
sampling_freq=44100
tone_freq=587
min_val=-2*np.pi
max_val=2*np.pi
生成时间轴和音频信号：
t=np.linspace(min_val, max_val, dutation*sampling_freq)
audio=np.sin(2*np.pi*tone_freq*t)
为信号增加一些噪声：
noise=0.4*np.np.random.rand(duration*sampling_freq)
audio+=noise
将这些值转换为16为整数，然后保存：
scaling_factoe=pow(2,15)-1
audio_normalized=audio/np.max(np.max(sudio))
audio_scaled=np.int16(audio_normalized*scaling_factor)
将信号写入输出文件：
write(output_file, sampling_freq, audio_scaled)
用前100个值画出该信号：
audio=audio[:100]
生成时间轴：
x_values=np.arange(0, len(audio), 1)/float(sampling_freq)
将时间轴的单位转换为秒：
x_values*=1000
画出信号：
plt.plot(x_values, audio, color='black')
plt.xlabel('Time (ms)')
plt.ylabel('Amplitude')
plt.title('Audio signal')
plt.show()
4. 合成音乐：
import json
import numpy as np
from scipy.io.wavfile import write
import matplotlib.pyplot as plt
定义一个函数，基于输入参数合成音调：
def synthesizer(freq, duration, amp=1.0, sample_freq=44100):
t=np.linspace(0, duration, duration*sampling_freq) # 创建时间轴
audio=amp*np.sin(2*np.pi*freq*t) # 构建音频信号
return audio.astype(np.int16)
定义main函数。json文件中包含音阶及频率。
if __name__=='__main__':
tone_map_file='tone_freq_map.json'
with open(tone_map_file,'r') as f:
tone_freq_map=json.loads(f.read())
tone imput_tone='G' # G调
duration=2
ampltude=10000
sampling_freq=44100
synthesized_tone=synthesized(tone_freq_map[input_tone], duration, amplitude, sampling_freq) # 生成音阶
write('output_tone.wave', sampling_freq, synthesized_tone) # 写入文件
tone_seq=[('D',0.3), ('G',0.6), ('C',0.5), ('A',0.3), ('Asharp',0.7)] # 音阶及持续时间
output=np.arrsy([])
for item in tone_seq:
input_tone=item[0]
duration=item[1]
synthesized_tone=synthesized(tone_freq_map[input_tome], duration, amplitude, s ampling_freq)
output=np.append(output, synthesized_tone, axis=0)
write('output_tone_seq.wav', samping_freq, output)
5. 提取频域特征：
在多数的现代语音识别系统中，会用到频域特征。将信号转换为频域后，还需要将其转换成有用的形式。梅尔频率倒谱系数MFCC（Mal Frequency Cepstrum Coefficient）首先计算信号的功率谱，然后用滤波器组和离散余弦变换的组合来提取特征。
需要安装python_speech_features包。
import numpy as np
import matplotlib.pyplot as plt
from scipy.io import wavfile
from features import mfcc,logfbank
读取文件：
sampling_freq,audio=wavfile.read('input_freq.wav')
提取MFCC和过滤器组特征：
mfcc_features=mfcc(audio, sampling_freq)
filterbank_features=logfbank(audio, sampling_freq)
打印参数，查看可生成多少个窗体：
print('\nMFCC:\nNumber of windows=', mfcc_features.shape[0])
print('Lenght of each feature=', mfcc_features.shape[1])
print('\nFilter bank:\nNumber of windows=', filterbank_features.shape[0])
print('Length of each feature=', filterbank_features.shape[1])
将MFCC特征可视化，转换矩阵时时域是水平的：
mfcc_features=mfcc_features.T
plt.matshow(mfcc_features)
plt.title('MFCC')
将滤波器组特征可视化，也需要转换矩阵，使得时域是水平的：
filterbank_features=filterbank_features.T
plt.matshow(filterbank_features)
plt.title('Filter bank')
plt.show()
6. 创建隐马尔可夫模型类：
隐马尔可夫模型HMMs（Hidden Markov Nodels）非常擅长建立时间序列数据模型，需要使用hmmlearn包。
创建类处理HMM相关过程：
class HMMTrainer(object):
def __init__(self, model_name='GaussianHMM', n_compoments=4, cov_type='diag', n_iter=1000):
self.model_name=model_name
self.n_components=n_components # 隐藏状态的个数
self.cov_type=cov_type # 转移矩阵的协方差类型
self.n_iter=n_tier # 迭代次数
self.models=[]
if self.model_name=='GaussianHMM': # 定义参数模型
self.model=hmm.GaussianHMM(n_compinents=self.n_components, covariance_type=self.cov_type, n_iter=felf.n_iter)
else:
raise TypeError('Invalid model type')
def train(self, X):
np.seterr(all='ignore')
self.models.append(delf.model.fit(X))
def get_score(self, input_data):
return self.model.score(input_data)
类中创建了处理隐马尔可夫模型的训练和预测。
7. 创建一个语音识别器：
语音文件数据库包含7个不同的单词，并且每个单词都有15个音频文件与之相关，需要为每一类构建一个隐马尔科夫模型。如果想识别新的输入文件中的单词，需要对该文件运行所有的模型，并找出最佳分数的结果。其中使用前面构建的隐马尔科夫类。
import os
omport argparse
import numpy as np
from scipy.io import wavfile
from hmmlearn import hmm
from features import mfcc
定义函数解析命令行中输入的参数：
def build_arg_parse():
parser=argparse.ArgumentParse(decription='Trains the HMM classifier')
parser.add_argument('--input-folder', dest='input_folder', required=True, help='Input filder containing the audio files insubfolders')
return parser
定义main函数，解析输入参数：
if __name__=='__main__':
args=build_arg_parser().parse_args()
input_filder=args.input_folder
hmm_models=[]
for dirname in os.listdir(input_folder): # 解析输入路径
subfolder=os.path.join(input_folder, dirname) # 提取子文件夹名称
if not os.path.isdir(subfolder):
continue
label=subfolder[subfolder.rfinf('/')+1:] # 提取标记，子文件夹名称即为类标记
X=np.array([])
y_words=[]
for filename in [x for x in os.listdir(subfolder) if x.endswith('.wav')][:-1]: # 迭代
filepath=os.path.join(subfolder, filename)
sampling_freq, audio=wavfile.read(filepath) # 读取每个音频文件
mfcc_feartures=mfcc(audio, sampling_freq) # 提取mfcc特征
if len(X)==0: # 将MFCC特征添加到X变量
X=mfcc_features
else:
X=np.append(X, mfcc_features, axis=0)
y_words.append(label) # 添加标记
hmm_trainer=HMMTrainer() # 训练HMM模型
hmm_trainer.train(X)
hmm_models.append((hmm_trainer, label)) # 保存HMM模型
hmm_trainer=None
# 测试文件
input_files=['data/pineapple/pineapple15.wav',
'data/orange/orange15.wav'
'data/apple/apple15.wav'
'data/kiwi/kiwi15.wav']
for input_file in input_files: # 解析输入文件
sampling_freq,audio=wavfile.read(input_file)
mfcc_features=mfcc(audio, sampling_freq) # 提取MFCC特征
max_score=None # 定义变量
output_label=None
for item in hmm_models: # 迭代所有hmm模型，并通过每个模型运行输入文件
hmm_model,label=item
score=hmm_model.get_score(mfcc_features) # 提取分数
if scire>max_score # 保存最大分数
max_score=score
output_label=label
print('\nTrue: ', input_file[input_file.find('/')+1:input_file.rfind('/')]
print('Predicted: ', output_label)
九、解剖时间序列和时序数据：
时间序列数据就是随时间的变化收集的测量序列数据，这些数据是根据预定义的变量并在固定的间隔时间采集的。时间序列数据最主要特征就是其顺序。
1. 将数据转换为时间序列格式：
需要使用pandas库来分析时间序列数据，需要先安装pandas库。
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
定义一个函数来读取输入文件，该文件将序列观察结果转换为时间序列数据。数据包含4列，前两列为年月，后两列为数据。
def convert_data_to_timeseries(input_file, column, verbose=False):
datanp.loadtxt(input_file,德力米特人，‘）’
# 提取第一列的起始时间，最后一列的终止时间
start_date=str(int(data[0,0]))+'-'+str(int(data[0,1]))
end_date=str(int(data[-1,0]+1))+'-'+str(int(data[-1,1]%12+1))
if verbose:
print('\nStart date=', start_date)
print('\nEnd date=', end_date)
# 创建以月为间隔的变量
dates=pd.date_range(start_date, end_date, freq='M')
# 将给定的列转换为时间序列数据
data_timeseries=pd.Series(data[:,column], index=dates)
# 打印最开始的10个元素
if verbose:
print('\nTime series data:\n', data_timeseries[:10]
return data_timeseries
定义主函数：
if __name__='__main__':
input_file='data_timeseries.txt'
# 加载文本文件第3列，并将其转换为时间序列数据
column_num=2
data_timeseries=convert_data_to_timeseries(input_file, column_num)
# pandas库提供了画图功能，画出数据序列
data_timeseries.plot()
plt.title('Input data')
plt.show()
2. 切分时间序列数据：
使用pandas切分时间序列数据，从时间序列数据的各个阶段获得相应信息。
column_num=2
data_timeseries=convert_data_to_timeseries(input_file, column_num)
start='2008' # 起始年份
end='2015' # 终止年份
plt.figure()
data_timeseries[start:end].plot() # 画出给定年份范围内的数据
plt.title('Data from '+start+' to '+end)
start='2007-2' # 画图起始年月
end='2007-11'
plt.figure()
data_timeseries[start:end].plot()
plt.title('Data from '+start+' to '+end)
plt.show()
3. 操作时间序列数据：
pandas库提供了各种操作时间序列数据的方式。
data1=convert_data_to_timeseries(input_file, 2)
data2=convert_data_to_timeseries(input_file, 3)
dataframe=pd.DataFrame({'first':data1, 'second':data2}) # 转化为pandas数据帧
dataframe['1952':'1955'].plot() # 画出给定年份范围内的数据
plt.title('Data overlapped on top of each other')
plt.figure()
difference=dataframe['1952':'1955']['first']-dataframe['1952':'1955']['second']
difference.plot() # 画出两组数据的不同
plt.title('Difference(first-second)')
# 设定不同阈值来过滤数据并绘图
dataframe[(dataframe['forst']>60)&(dataframe['second']N20)].plot()
plt.title('first>60 and second<20')
plt.show()
4. 从时间序列数据中提取统计数字：
从时间序列数据中提取出有用的统计信息。
data1=convert_data_to_timeseries(input_file, 2)
data2=convert_data_to_timeseries(input_file, 3)
dataframe=pd.DataFrame({'first':data1, 'second':data2}) # 转化为pandas数据帧
print('\nMaximum:\n', dataframe.max()) # 打印最大值
print('\nMinimum:\n', dataframe.min()) # 打印最小值
print('\nMean:\n', dataframe.mean()) # 打印数据均值
print('\nMean row-wise:\n', dataframe.mean(1)[:10]) # 打印每行均值
# 滑动均值是指计算一个窗口内的信号均值，并不断地移动时间窗，窗口大小24
pd.rolling_mean(dataframe, window=24).plot()
print('\nCorrelation coefficients:\n', dataframe.corr()) # 打印相关性系数
# 画出滑动相关性，窗口大小60
plt.figure()
pd.rolling_corr(dataframe['first'], dataframe['second'], window=60).plot()
plt.show()
5. 针对序列创建隐马尔科夫模型：
隐马尔科夫模型HMMs是处理序列数据的强大方法，广泛应用于金融、语音分析、天气预测、单词序列等领域。要使用hmmlearn程序包。
import datetime
import numpy as np
import matplotlib.pyplot as plt
from hmmlearn.hmm import GaussianHMM
from convert_to_timeseries import convert_data_to_timeseries
使用的data文件bak带逗号分隔符的行，每行包括3个值，年份、月份、浮点数据。加载到NumPy数组。
input_file='data_hmm.txt'
data=np.loadtxt(input_file, delimiter=',')
如果数据有多列，可将数据按照列的方向堆叠起来用于分析。
X=np.column_stack([data[:,2]])
用4个成分创建并训练HMM。成分的个数是一个需要进行选择的超参数。这里选择4个成分，意味着用4个潜在状态生成数据。
print('\nTraining HMM......')
num_components=4
model=GaussianHMM(n_components=num_components, covariance_type='diag', n_iter=1000)
model.fit(X)
运行预测器以获得隐藏状态：
hidden_states=model.predict(X)
计算这些隐藏状态的均值和方差：
for i in range(model.n_components):
print('\nHidden state', i+1)
print('Mean=', round(model.means_[i][0], 3)
print('Variance=', round(np.diag(model.covars_[i])[0], 3)
HMM是一个生成模型，这里生成1000个示例数据并将其画出：
num_samples=1000
samples,_=model.sample(num_samples)
plt.plot(np.arange(num_samples), samples[:,0], c='black')
plt.title('Number of compinents='+str(num_components))
plt,show()
6. 针对序列文本数据创建条件随机场：
条件随机场CRFs（Conditional Random Field）是一个概率模型，用于分析结构化数据，常用于标记和分段序列数据。这是一个判断模型，用于分析序列、股票、语音、单词等。在模型中，给定一个带标签的观察序列，对这个序列定义一个条件随机分布。使用pystruct库。
import os
inport argparse
import cPickle as pickle
import numpy as np
import matplotlib.pyplot as plt
from pystruct.datasets import load_letters
from pystruct.models import ChainCRF
from pystruct.learners import FrankWolfeSSVM
定义一个参数解析器，并将C值作为输入参数。C是一个超参数，用于控制模型：
def build_arg_parser():
parser=argparse.ArgumentParser(description='Trains the CRF classifier')
parser.add_argument('--c-value', dest='c_value', required=False, type=float, default=1.0, help='The C value that will be used for training')
return parser
定义一个处理所有与CRF相关处理的类：
class CRFTrainer(object):
def __init__(self, c_value, classifier_name='ChainCRF'):
self.c_value=c_values
self.classifier_name=classifier_name
if self.classifier_name=='ChainCRF': # 错误检查
model=ChainCRF() # 用到分类器，下面使用SVM
self.clf-FrankWolfeSSVM(model=model, C=self.c_value, max_iter=50)
else:
raise TypeError('Invalid classifier type')
# 加载字母数据集，其中包括分割的字母以及与其相关的特征向量
def load_data(self):
letters=load_letters()
X,y,folds=letters['data'], letters['labels'], letters['folds'] # 加载数据和标签到相应变量
X,y=np.array(X), np.array(y)
return X,y,folds
# 定义一个训练方法：
def train(self, X_train, y_train):
self.clf.fit(X_train, y_train)
# 定义一个方法来评价模型的性能
def evaluate(self, X_test, y_test):
return self.clf.score(X_test, y_test)
# 定义一个方法，将新数据进行分类：
def classifier(self, input_data):
return self.clf.predict(input_data)[0]
字母的编号在数组中被索引，为了可读，定义一个函数将数字转换为字母。
def dicoder(arr):
alphabets='abcdefghijklmnopqrstuvwxyz'
output=''
for i in arr:
output+=alphabets[i]
return output
定义主函数并解析输入参数：
if __name__=='__main__':
args=build_arg_parse().parse_args()
c_value=args.c_value
# 用类和C值初始化变量
crf=CRFTrainer(c_value)
# 加载字母数据
X, y, folds=crf.load_data()
# 将数据分成训练数据集和测试数据集
X_train, X_test=X[folds==1], X[folds!=1]
y_train, y_test=y[folds==1], y[folds!=1]
# 训练CRF模型
crf.train(X_train, y_train)
# 评价CRF模型的性能
score=crf.evaluate(X_test, y_test)
print('\nAccuracy score=', str(round(score*100, 2))+'%'
# 输入一个随机测试向量，并用这个模型预测输出
print('\nTrue label=', decoder(y_test[0])
predicted_output=crf.classify([X_test[0]])
print('Predicted output=', decoder(predicted_output)
7. 用隐马尔科夫模型分析股票市场数据：
股票市场数据是典型的时间序列数据。
import datetime
import numpy as np
import matplotlib.pyplot as plt
import matplotlib.finance import quotes_historical_yahoo_ochl
from hmmlearn.hmm import GaussianHMM
从雅虎财经获取股票报价：
quotes=quotes_historical_yahoo_ichl('INTC', datetime.date(1994,4,5), datetime.date(2015,7,3))
每个报价包含6个值，提取需要的数据：
dates=np.array([quote[0] for quote in quotes], dtype=np.int)
closing_values=np.array([quote[2] for quote in quotes])
volume_of_shares=np.array([quote[5] for quote in quotes])[1:]
计算每天收盘价的变化率，用这个变化率作为一个特征：
diff_percentage=100.0*np.diff(closing_values)/closing_values[:-1]
dates=dates[1:]
将数组进行列堆叠以用作训练：
X=np.column_stack([diff_percentage, volume_of_shares])
用5个成分训练隐马尔科夫模型：
model=GaussianHMM(n_components=5, coveriance_type='diag', n_iter=1000)
model.fit(X)
生成500个示例数据用于训练隐马尔科夫模型并画出：
num_samples=500
samples,_=model.sample(num_samples)
plt.plot(np.arange(num_samples), samples[:,0], c='black')
plt.show()
十、图像内容分析：
计算机视觉是一个研究如何处理、分析和理解视觉数据内容的领域，包括目标识别、形状分析、姿态估计、3D建模、视觉搜索等。
计算机视觉包括多个级别的分析，在低级视觉分析领域，可以进行像素处理，如边检测、形态处理和光流；在中高级分析领域，可以处理事务，如物体识别、建模、运动分析等。高层次分析往往依赖低层次分析的输出结果。
OpenCV是世界上最受欢迎的计算机视觉库，已成业界的事实标准，使用此库来分析图像。
1. 用OpenCV-Python操作图像：
import sys
import cv2
import numpy as np
指定输入图像为文件的第1个参数，并使用图像读取函数来读取参数：
input_file=sys.argv[1]
img=cv2.imread(input_file)
显示输入图像：
cv2.imshow('Original', img)
剪裁图像，提取图像的高度和宽度，然后指定边界：
h,w=img.shape[:2]
start_row,end_row=int(0.21*h), int(0.73*h)
start_col,end_col=int(0.37*w), int(0.92*w)
用NumPy切分方式剪裁图像，并将其展示出来：
img_cropped=img[start_row:end_row, start_col:end_col]
cv2.inshow('Cropped', img_cropped)
将图像大小调整为其原始大小的1.3倍，并将其展示出来：
scaling_factor=1.3
img_scaled=cv2.resize(img, None, fx=scaling_factor, fy=Scaling_factor, interpolation=cv2.INTER_LINEAR)
cv2.imshow('Uniform resizing', img_scaled)
如果只希望在一个维度上进行调整：
img_scaled=cv2.resize(img, (250,400), interpolation=cv2.INTER_AREA)
cv2.imshow('Skewed resizing', img_scaled)
将图像保存到输出文件：
output_file=input_file[:-4]+'_cropped.jpg'
cv2.imwrite(output_file, img_cropped)
cv2.waitKey()
waitKey函数保持显示图像，直到按下键盘上的任一按键。
2. 检测边：
边缘检测是计算机视觉中最常用到的技术之一。
import sys
import cv2
import numpy as np
加载输入图像，并转换为灰度图像：
input_file=sys.argv[1]
img=cv2.imread(input_file, cv2.IMREAD_GRAYSCALE)
提取输入图像的高度和宽度：
h,w=img.shape
Sobel滤波器是一种边缘检测器，使用3x3内核来检测水平边和垂直边：
sobel_horizontal=cv2.Sobel(img, cv2.CV_64F, 1, 0, ksize=5) # Sobel水平检测
sobel_vertical=cv2.Sobel(img, cv2.CV_64F, 0, 1, ksize=5) # Sobel垂直检测
拉普拉斯边缘检测器可以检测两个方向上的边：
laplacian=cv2.Laplacian(img, cv2.CV_64F)
Canny边缘检测器在解决噪声问题方面由于拉普拉斯检测器，这是一个分阶段的处理过程：
canny=cv2.Canny(img, 50, 240)
显示所有的输出图像：
cv2.imshow('Original', img)
cv2.imshow('Sobel horizontal', sobel_horizontal)
cv2.imshow('Sobel vertical', sobel_vertical)
cv2.imshow('Laplacian', laplacian)
cv2.imshow('Canny', canny)
cv2.waitKey()
3. 直方图均衡化：
指修改图像的像素以增强图像的对比度强度。
import sys
import cv2
import numpy as np
加载输入图像：
input_file=sys.argv[1]
img=cv2.imread(input_file)
转化为灰度图像：
img_gray=cv2.cvtColor(img, cv2.COLOR_BGR2GRAY)
cv2.inshow('Input grayscale image', img_gray)
均衡灰度图像的直方图，并显示出来：
img_gray_histeq=cv2.equalizeHist(img_gray)
cv2.inshow('Histogram equalized-grayscale', img_gray_histeq)
直方图均衡化只适用于亮度通道。如果对彩色图像，需要先转化为YUV色彩空间，均衡Y通道，然后再将其转换回RGB彩色图像。
img_yuv=cv2.cvtColor(img, cv2.COLOR_BGR2YUV)
img_yuv[:,:,0]=cv2.equalizeHist(img_yuv[:,:,0]) # 均衡Y通道
img_histeq=cv2.cvtColor(img_yuv, cv2.COLOR_YUV2BGR)
显示输入和输出图像：
cv2.imshow('Input color image', img)
cv2.inshow('Histogram equalized-color', img_histeq)
cv2.waitKey()
4. 检测棱角：
棱角检测帮助识别图像中突出的点，是一种早期的特征提取技术之一。
import sys
import cv2
import numpy as np
加载输入图像：
input_file=sys.argv[1]
img=cv2.imread(input_file)
cv2.imshow('Input image', img)
将图像转为灰度，并将其强制转换为浮点值：
img_gray=cv2.cvtColor(img, cv2.COLOR_BGR2GRAY)
img_gray=np.float32(img_gray)
对灰度图像运行哈里斯角检测器（Harris corner detector）函数：
img_harris=cv2.cornerHarris(img_gray, 7, 5, 0.04)
为了标记棱角，需要放大图像：
img_harris=cv2.dilate(img_harris, None)
定义显示重要点个数的阈值：
img[img_harris>0.01*img_harris.max()]=[0, 0, 0]
显示输出图像：
cv2.inshow('Harris Corners', img)
cv2.waitKey()
5. 检测SIFT特征点：
尺度不变特征变换SIFT（Scale Invariant Feature Transform）是计算机视觉领域最常用的特征之一，其在大小、方向、对比度方面都有较强健壮性。
import sys
import cv2
import numpy as np
加载输入图像：
input_file=sys.argv[1]
img=cv2.imread(input_file)
将图像转换为灰度：
img_gray=cv2.cvtColor(img, cv2.COLOR_BGR2GRAY)
初始化SIFT检测器对象并提取关键点：
sift=cv2.xfeatures2d.SIFT_create()
keypoints=dift.detect(img_gray, None)
在输入图像上画出关键点：
img_sift=np.copy(img)
cv2.drawKeypoints(img, keypoints, img_sift,
flags=cv2.DRAW_MATCHES_FLAGS_DRAW_RICH_KEYPOINTS)
显示输入和输出图像：
cv2.imshow('Input image', img)
sv2.imshow('SIFT features', img_sift)
cv2.waitKey()
6. 创建Star特征检测器：
当创建目标识别系统时，在用SIFT检测特征之前，可能需要用到一个不同的特征检测器，能通过灵活地层叠不同的模块来获得最佳的性能。
import sys
import cv2
import numpy as np
定义一个类，用于处理与Star特征检测相关的函数：
class StarFeatureDetector(object):
def __init__(self):
self.detector=cv2.xfeatures2d.StarDetector_create()
# 定义一个对输入图像运行检测器的函数
def detect(self, img):
return self.detector.detect(img)
在main主函数中加载输入图像：
if __name__=='__main__':
imput_file=sys.argv[1]
img=cv2.imread(input_file)
# 将图像转为灰度
img_gray=cv2.cvtColor(input_img, cv2.COLOR_BGR2GRAY)
# 用Star特征检测器检测出特征
keypoints=StarFeatureDetector().detect(input_img)
# 画出输入图像的关键点
cv2.drawKeypoints(input_img, keypoints, input_img,
flags=cv2.DRAW_MATCHES_FLAGS_DRAW_RICH_KEYPOINTS)
# 显示输出图像：
cv2.imshow('Star features', input_img)
cv2.waitKey()
7. 利用视觉码本和向量量化创建特征：
为了创建一个目标识别系统，需要从每张图像中提取特征向量。每张图像需要有一个识别标志，以用于匹配。用一个称视觉码本的概念来创建图像识别标志。在训练数据集中，这个码本基本上是一个字典，用于提出关于图像的描述。用向量量化方法将很多特征点进行聚类并得出中心点，这些中心点将作为视觉码本的元素。
为了创建一个健壮的目标识别系统，需要数万幅图像，有一个非常著名的数据集Caltech256，包括256类图像，每一类包含上千幅示例图像。
定义一个提取特征的类：
class FeatureBuilder(object):
定义一个从输入图像提取特征的方法，用Star检测器获得关键点，然后用SIFT提取这些位置的描述信息：
def extract_features(self, img):
keypoints=StarFeatureDetector().detect(img)
ketpoints,feature_vectors=compute_sift_features(img, keypoints)
return feature_vectors
从描述信息中提取出中心点：
def get_codewords(self, input_map, scaling_size, max_samples=12):
keypoints_all=[]
count=0
cur_class=''
每幅图像都会生成大量的描述信息，代码中仅用一小部分图像。
for item in input_map:
if count>=max_samples:
if cur_class!=item['object_class']:
count=0
else:
continue
count+=1
将进程打印出来：
if count==max_samples:
print('Built centroids for ', item['object_class']
# 提取当前标签
cur_class=item['object_class']
# 读取图像并调整其大小
img=cv2.imread(item['image_path'])
img=resize_image(img, scaling_size)
# 设置维度数为128并提取特征
num_dims=128
feature_vectors=self.extract_features(img)
keypoints_all.extend(feature_vectors)
# 用向量量化特征点
kmeans, centroids=BagOfWords().cluster(keypoints_all)
return kmeans, centroids
定义一个类处理词袋模型和向量量化：
class BagOfWords(object):
def __init__(self, num_clusters=32):
self.num_dims=128
self.num_clusters=num_clusters
self.num_retries=10
定义一个方法来量化数据点，用k-means聚类来实现：
def cluster(self, datapoints):
kmeans=KMeans(self.num_clusters, n_init=max(self.num_retries, 1), max_iter=10, tol=1.0)
# 提取中心点
res=knean.fit(datapoints)
centroids=res.cluster_centers_
return kmeans, centroids
定义一个方法来归一化数据：
def normalize(self, input_data):
sum_input=np.sum(input_data)
if sum_input>0:
return input_data/sum_input
else:
return input_data
定义一个方法来获得特征向量：
def construct_feature(self, img, kmeans, centroids):
keypoints=StarFeatureDetector().detect(img)
keypoints, feature_vectors=compute_sift_features(img, keypoints)
labels=kmeans.predict(feature_vectors)
feature_vector=np.zeros(self.num_clusters)
创建一个直方图并将其归一化：
for i,item in emumerate(feature_vectors):
feature_vector[labels[i]]+=1
feature_vector_img=np.reshape(feature_vector, ((1, feature_vector.shape[0])))
return self.normalize(feature_vector_img)
定义一个方法来提取SIFT特征：
def compute_sift_features(img, keypoints):
if img is None:
raise TypeError('Invalid input image')
img_gray=cv2.cvtColor(img, cv2.COLOR_BGR2GRAY)
keypoints, descriptors=cv2.xfeatures2d.SIFT_create().compute(img_gray, keypoints)
return keypoints, descriptors
8. 用极端随机森林训练图像分类器：
极端随机森林ERF（Extremely Random Forests）在机器学习中很流行，具有较快的速度和比较精确的准确度。基于图像的特征构建一组决策树，并通过训练这个森林实现正确决策。
import argparse
import cPickle as pickle
import numpy as np
from sklearn.ensemble import ExtraTreesClassifier
from sklearn import preprocessing
定义一个参数解析器：
def build_arg_parser():
parser=argparse.ArgumentParser(description='Trains the classifier')
parser.add_argument('--feature-map-file', dest='feature_map_file', required=True, help='Input pickle file containing the feature map')
parser.add_argument('-model-file', dest='model_file', required=False, help='Output file where the trained model will be stored')
return parser
定义一个类来处理ERF训练，用到一个标签编码器来对训练标签进行编码：
class ERFTrainer(object):
def __init__(self, X, label_words):
self.le=preprocessing.LabelEncoder()
self.clf=ExtraTreesClassifier(n_estimators=100, max_depth=16, random_state=0)
# 对标签编码并训练分类器
y=self.encode_labels(label_words)
self.clf.fit(np.asarray(X), y)
定义一个函数，用于对标签进行编码：
def classify(self, X):
label_nums=self.clf.predict(np.asarray(X))
label_words=self.le.inverse_transform([int(x) for x in label_nums])
return label_words
定义main函数并解析输入参数：
if __name__=='__main__':
args=build_arg_parser().parse_args()
feature_map_file=args.feature_map_file
model_file=args.model_file
# 加载特征地图
with open(feature_map_file, 'r') as f:
feature_map=pickle.load(f)
# 提取特征向量
label_words=[x['object_class'] for x in feature_map]
dim_size=feature_map[0]['feature_vector'].shape[1]
X=[np.reshape(x['feature_vector'], (dim_size,)) for x in feature_map]
# 基于训练数据训练ERF
erf=ERFTrainer(X, label_words)
# 保存训练的ERF模型
if rgs.model_file:
with open(args.model_file, 'w') as f:
pickle.dump(erf, f)
结果产生一个erf.pkl文件。
9. 创建一个对象识别器：
训练好一个ERF模型后，接着创建一个目标识别器，可用来识别未知图像的内容。
import argparse
import cPickle as pickle
import cv2
import numpy as np
import build_features as bf
from trainer import ERFTrainer
定义一个参数解析器：
def build_arg_parser():
parser=argparse.ArgumentParser(description='Extracts features from each line and classifier the data')
parser.add_argument('--input-image', dest='input_image', required=True, help='Input image to be classified')
parser.add_argument('-model-file', dest='model_file', required=True, help='Input file containing the trained model')
parser.add_argument('codebook-file', dest='codebook_file', required=True, help='Input file containing the codebook')
return parser
定义一个类来处理图像标签提取函数：
class ImageTagExtractor(object):
def __init__(self, model_file, cookbook_file):
with open(model_file, 'r') as f:
self.erf=pickle.load(f)
with open(cookbook_file, 'r') as f:
self.kmeans, self.centroids=pickle.load(f)
定义一个函数，用于使用训练好的ERF模型来预测输出：
def predict(self, img, scaling_size):
img=bf.resize_image(img, scaling_size)
feature_vector=bf.BagOfWords().construct_feature(img, self.kmeans, self.centroids)
image_tag=self.erf.classify(feature_vector)[0]
return image_tag
定义main函数，加载输入图像：
if __name__=='__main__':
args=build_arg_parser().parse_args()
model_file=args.model_file
codebook_file=args.cookbook_file
input_image=cv2.imread(args.input_image)
# 合理调整图像大小
scaling_size=200
# 打印输出结果
print('\nOutput:', OmageTagExtractor[model_file, codebook_file).predict(input_image, scaling_size)
十一、人脸识别：
人脸识别是指在给定图像中识别某个人；人脸检测是从给定图像中定位人脸位置。在一个典型的生物人脸识别系统中，需要在识别脸部之前确定脸部的位置。
1. 从网络摄像头采集和处理视频信息：
import cv2
OpenCV提供了一个视频采集对象，可以用来从摄像头采集图像，0参数指定网络摄像头ID：
cap=cv2.VideoCapture(0)
定义网络摄像头采集的帧的比例系数：
scaling_factor=0.5
启动一个无限循环来采集帧，直到按下Esc键：
while True:
ret,frame=cap.read()
# 调整帧的大小，可选：
frame=cv2.resize(frame, None, fx=scaling_factor, fy=scaling_factor, interpolation=cv2.INTER_AREA)
# 显示帧
cv2.imshow('Webcam', frame)
# 等待1ms，然后采集下一帧
c=cv2.waitKey(1)
if c==27:
break
# 释放视频采集对象
cap.release()
# 在结束代码之前关闭所有获得窗体
cv2.destoryAllWondows()
2. 用Haar级联创建一个人脸识别器：
用Haar级联做人脸检测。Haar级联通过在多个尺度上从图像中提取大量的简单特征来实现。简单特征主要指边缘、线、矩形特征等，这些特征都非常容易计算，然后通过创建一系列简单的分类器来做训练，使用自适应增强技术可以使得这个过程更健壮。
import cv2
import numpy as np
加载人脸检测级联文件，这是可以用作检测器的训练模型：
face_cascade=cv2.CascadeClassifier('cascade_files/haarcascade_frontalface_alt.xml')
确定级联文件是否正确地加载：
if face_cascade.empty():
raise IOError('Unable to load the face cascase classifier xml file')
生成一个视频采集对象：
cap=cv2.VideoCapture(0)
定义图像向下采用的比例系数：
scaling_factor=0.5
循环采集直到按下Esc键：
while True:
ret,frame=cap.read()
# 调整帧的大小
frame=cv2.resize(frame, None, fx=scaling_factor, fy=scaling_factor,
interpolation=cv2.INTER_AREA)
# 将图像转换为灰度，使用灰度图像进行人脸检测
gray=cv2.cvColor(frame, cv2.COLOR_BGR2GRAY)
'''在灰度图像上进行人脸检测，参数1.3是指每个阶段的乘积系数，参数5是指每个候选矩形应该拥有的最小近邻数量。候选矩形是指人脸可能被检测到的候选区域'''
face_rects=face_cascade.detectMultiScale(gray, 1.3, 5)
# 对于每个检测到的人脸区域，在其周围画出矩形
for (x,y,w,h) in face_rects:
cv2.rectangle(frame, (x,y), (x+w,y+h), (0,255,0), 3)
# 展示输出图像
cv2.imshow('Face Detector', frame)
# 下下一次迭代前等待1ms，如果用户按下Esc键，退出循环
c=cv2.waitKey(1)
if c==27:
break
在结束代码前，释放并销毁所有对象：
cap.release()
cv2.destroyAllWindows()
3. 创建一个眼睛和鼻子检测器：
Haar级联方法可以被扩展应用于各种对象的检测。
import cv2
import numpy as np
加载人脸、眼睛和鼻子级联文件：
face_cascade=cv2.CascadeClassifier('cascade_files/haarcascade_frontalface_alt.xml')
eye_cascade=cv2.CascadeClassifier('cascade_files/haarcascade_eye.xml')
nose_cascade=cv2.CascadeClassifier('cascade_files/haarcascade_nose.xml')
确定文件释放正确加载：
if face_cascade.empty():
raise IOError('Unable to load the face cascase classifier xml file')
if eye_cascade.empty():
raise IOError('Unable to load the eye cascase classifier xml file')
if nose_cascade.empty():
raise IOError('Unable to load the nose cascase classifier xml file')
初始化视频采集对象：
cap=cv2.VideoCapture(0)
定义比例系数：
scaling_factor=0.5
重复循环直至用户按下Esc键：
while True:
ret,frame=cap.read()
# 调整帧的大小
frame=cv2.resize(frame, None, fx=scaling_factor, fy=Scaling_factor, interpolation=cv2.INTER_AREA)
# 将图像转换为灰度图像：
gray=cv2.cvtColor(frame, cv2.COLOR_BGR2GRAY)
# 在灰度图像上进行人脸检测
faces=face_cascade.detectMultiScale(gray, 1.3, 5)
# 在脸部区域运行眼睛和鼻子检测
for (x,y,w,h) in faces:
# 提取人脸ROI信息
roi_gray=gray[y:y+h, x:x+w]
roi_color=frame[y:y+h, x:x+w]
# 运行眼睛、鼻子检测
sys_rects=eye_cascade.detectMultiScale(roi_gray)
nose_rects=nose_cascade.detectMultiScale(roi_gray, 1.3, 5)
# 在眼睛周围画绿色的窗
for (x_eye, y_eye, w_eye, h_eye) in eye_rects:
center=(int(x_eye+0.5*w_eye), int(y_eye+0.5*h_eye))
radius=int(0.3*(w_eye+h_eye))
color=(0,255,0)
thickness=3
cv2.circle(roi_color, center, radius, color, thickness)
# 在鼻子周围画矩形
for (x_nose, y_nose, w_nose, h_nose) in nose_rects:
cv2.rectangle(roi_color, (x_nose,y_nose), (x_nose+w_nose, y_nose+h_nose), (0,255,0), 3)
break
# 展示图像
cv2.imshow('Eye and nose detector', frame)
# 在下一次迭代前等待1ms，如果用户按下Esc键，就跳出循环
c=cv2.waitKey(1)
if c==27:
break
结束代码前释放并销毁所有对象：
cap.release()
cv2.destroyallWindows()
4. 主成分分析：
主成分分析PCA（Principal Components Analysis）是一个降低维度的技术，常用于计算机视觉和机器学习中。当需要处理很大的特征维度时，训练机器学习系统变得异常困难，因此需要在训练之前降低数据的维度，但降低维度时并不想损失数据中的重要信息。PCA识别数据中的重要成分，并将其按照重要程度排序。
import numpy as np
from sklearn import decomposition
为输入数据定义5个维度，前两个维度是独立的，但后3个维度将依赖于前两个维度，即可以去掉后3个维度，因为并没有提供任何新信息：
x1=np.random.normal(size=250)
x2=np.random.normal(size=250)
x3=2*x1+3*x2
x4=4*x1-x2
x5=x3+2*x4
创建一个带这些特征的数据集：
X=np.c_[x1, x3, x2, x5, x4]
创建一个PCA对象：
pca=decomposition.PCA()
对输入数据做主成分分析：
pca.fit(X)
打印维度的方差：
variances=pca.explained_variance_
print('\nVariances in decreasing order:\n', variances)
如果一个特定的维度是有用的，那么它应该有一个有意义的方差值。设置一个阈值并确定重要的维度：
thresh_variance=0.8
num_useful_dims=len(np.where(variances>thresh_variance)[0])
print('\nNumber of userful dimensions:', num_useful_dims)
PCA识别只有两个维度在这个数据集中很重要：
pca.n_components=num_useful_dims
将数据集从5维转换为二维：
X-new=pca.fit_transform(X)
print('\nShape before:', X.shape)
print('Shape after:', X_new.shape)
5. 做核主成分分析：
PCA是以线性方式工作的，如果数据集不是以线性方式组织的，PCA就不能实现其功能。而核主成分分析则可以很好地解决这个问题。
import numpy as np
import matplotlib.pyplot as plt
from sklearn.decomposition import PCA, KernelPCA
from sklearn.datasets import make_circles
定义随机数发生器的种子值，以产生用于分析的数据示例：
np.random.seed(7)
生成以同心圆分布的数据：
X,y=make_circles(n_samples=500, factor=0.2, noise=0.04)
对这组数据做主成分分析：
pca=PCA()
X_pca=pca.fit_transform(X)
对输入数据做核主成分分析：
kernel_pca=KernelPCA(kernel='rbf', fit_inverse_transform=True, gamma=10)
X_kernel_pca=kernel_pca.fit_transform(X)
X_inverse=kernel_pca.inverse_transform(X_kernel_pca)
画出原始输入数据：
class_0=np.where(y==0)
class_1=np.where(y==1)
plt.figure()
plt.title('Original data')
plt.plot(X[class_0,0],X[class_0,1],'ko',mfc='none')
plt.plot(X[class_1,0],X[class_1,1],'kx')
plt.xlabel('1st dimensiton')
plt.ylabel('2nd dimensiton')
画出主成分分析后的数据：
plt.figure()
plt.plot(X_pca[class_0,0],X_pca[class_0,1],'ko',mfc='none')
plt.plot(X_pca[class_1,0],X_pca[class_1,1],'kx')
plt.title('Data transformed using PCA')
plt.xlabel('1st principal component')
plt.ylabel('2nd principal component')
画出核主成分分析后的数据：
plt.figure()
plt.plot(X_kernel_pca[class_0,0],X_kernel_pca[class_0,1],'ko',mfc='none')
plt.plot(X_kernel_pca[class_1,0],X_kernel_pca[class_1,1],'kx')
plt.title('Data transformed using Kernel PCA')
plt.xlabel('1st principal component')
plt.ylabel('2nd principal component')
用核方法将数据转换回原始空间，查看是否存在可逆关系：
plt.figure()
plt.plot(X_inverse[class_0,0],X_inverse[class_0,1],'ko',mfc='none')
plt.plot(X_inverse[class_1,0],X_inverse[class_1,1],'kx')
plt.title('Inverse transform')
plt.xlabel('1st dimensiton')
plt.ylabel('2nd dimensiton')
plt.show()
6. 盲源分离：
盲源分析是指将信号从混合体中分离出来。假设一组不同的信号发生器生成了不同的信号，而一个公共接收机接收到了所有这些信号。利用这些信号的性质，将这些信号从混合体中分离出来。使用独立成分分析ICA（Independent Components Analysis）算法来实现。
import numpy as np
import matplotlib.pyplot as plt
from scipy import signal
from sklearn.decomposition import PCA,FastICA
加载txt文件中的数据：
input_file='mixture_of_signals.txt'
X=np.loadtxt(input_file)
创建ICA对象：
ica=FastICA(n_components=4)
基于ICA重构信号：
signals_ica=ica.fit_transform(X)
提取混合矩阵：
mixing_mat=ica.mixing_
执行PCA做对比：
pca=PCA(n_components=4)
signals_pca=pca.fit_transform(X)
定义信号列表来将其画出：
models=[X, signals_ica, signals_pca]
指定颜色：
colors=['blue', 'red', 'black', 'green']
画出输入信号：
plt.figure()
plt.title('Input signal (mixture)')
for i, (sig,color) in enumerate(zip(X.T, colors), 1):
plt.plot(sig, color=color)
画出利用ICA分离的信号：
plt.figure()
plt.title('ICA serarated signals')
plt.subplots_adjust(left=0.1, bottom=0.05, right=0.94, top=0.94, wspace=0.25, hspace=0.45)
用不同的颜色画出子图：
for i,(sig,color) in enumerate(zip(signals_ica.T, colors), 1):
plt.subplot(4,1,i)
plt.title('Signal '+str(i))
plt.plot(sig, color=color)
画出用PCA分离的信号：
plt.figure()
plt.title('PCA serarated signals')
plt.subplots_adjust(left=0.1, bottom=0.05, right=0.94, top=0.94, wspace=0.25, hspace=0.45)
用不同的颜色画出各子图：
for i,(sig,color) in enumerate(zip(signals_pca.T, colors), 1):
plt.subplot(4,1,i)
plt.title('Signal '+str(i))
plt.plot(sig, color=color)
plt.show()
7. 用局部二值模式直方图创建一个人脸识别器：
import os
import cv2
import numpy as np
from sklearn inport preprocessing
定义一个类来处理与类标签编码相关的所有任务：
class LabelEncoder(object):
定义一个方法来为这些标签编码。在输入训练数据中，标签用单词表示，需要转为数字来训练系统。
def encode_labels(self, label_words):
self.le=preprocessing.LabelEncoder()
self.le.fit(label_words)
定义一个将单词转换成数字的方法：
def word_to_num(self, label_word):
return int(self.le.transform([label_word])[0])
定义一个方法，用于将数字转换会其原始单词：
def num_to_word(self, label_num):
return self.le.inverse_transform([label_num])[0]
定义一个方法，用于从输入文件夹中提取图像和标签：
def get_images_and_labels(input_path):
label_words=[]
# 对输入文件夹做递归迭代，提取所有图像的路径
for root,dirs,files in os.walk(input_path):
for filename in (x for x in files if x.endswitch('.jpg')):
filepath=os.path.join(root, filename)
label_words.append(filepath.split('/')[-2])
# 初始化变量
images=[]
le=LabelEncoder()
le.encode_labels(label_words)
labels=[]
# 为训练解析输入目录
for root,dirs,files in os.walk(input_path):
for filename in (x for x in files if x.endswith('.jpg')):
filepath=os.path.join(root,filename)
# 将当前图像读取成灰度格式
image=cv2.imread(filepath,0)
# 从文件夹路径中提取标签
name=filepath.split('/')[-2]
# 对该图像做人脸检测
faces=faceCascade.detectMultiScale(image, 1.1, 2, minSize=(100,100))
# 提取ROI属性值，并将这些值和标签编码器返回
for (x,y,w,h) in faces:
images.append(image[y:y+h, x:x+w])
labels.append(le.word_to_num(name))
return images,labels,le
定义main函数，并定义人脸级联文件的路径：
if __name__=='__main__':
cascade_path='cascade_files/haarcascade_frontalface_alt.xml'
path_train='faces_dataset/train'
path_test='faces_dataset/test'
# 加载人脸级联文件
faceCascade=cv2.CascadeClassifier(cascade_path)
# 生成局部二值模式直方图人脸识别器对象
recognizer=cv2.face.createLBPHFaceRecognizer()
# 为输入路径提取图像、标签和标签编码器
images,labels,le=get_images_and_labels(path_train)
# 用提取的数据训练人脸识别器
recognizer.train(images, np.array(labels))
# 用未知数据测试人脸识别器
stop_flag=False
for root,dirs,files in os.walk(path_test):
for filename in (x for x in files if x.endwith('.jpg')):
filepath=os.path.join(root, filename)
# 加载图像
predict_image=cv2.inread(filepath, 0)
# 用人脸检测器确定人脸的位置
faces=faceCascade.detectMultiScale(predict_image, 1.1, 2, minSize=(100,100))
# 对于每一个人脸ROI，运行人脸识别器
for (x, y, w, h) in faces:
predicted-index, cinf=recognizer.predict(predict_image[y:y+h, x:x+w])
# 将标签转换为单词
predicted_person=le.num_to_word(predicted_index)
# 在输出图像中叠加文字，并将其展示
cv2.putText(predict_image, 'Prediction:'+predicted_person, (10,60). cv2.FONT_HERSHEY_SIMPLEX, 2, (255,255,255), 6)
cv2.imshow('Recignizing face', predict_image)
# 检查用户是否按下Esc键，如有跳出循环
c=cv2.waitKey(0)
if c==27:
stop_flag=True
break
if stop_flag:
break
十二、深度神经网络：
神经网络就是模仿人类大脑激发学习过程的框架。神经网络处理的是数字，如果要处理包含图像、文字、传感器等的任务，就必须将其转换成数值形式，然后将其输入到一个神经网络。可以用神经网络做分类、聚类、生成以及其他相关的任务。
神经网络由一层层神经元组成，这些神经元模拟人类大脑中的生物神经元，每一层都是一组独立的神经元，这些神经元与相邻层的神经元相连。输入层对应提供的输入数据，而输出层包括了期望的输出结果，输入层与输出层之间的层统称为隐藏层。
为了使神经网络完成相应的任务，需要提供带标签的训练数据。神经网络将通过优化成本函数来训练自己，不停地迭代，直到错误率下降到一定的阈值为止。
深度神经网络是由多个隐藏层组成的神经网络，深度学习用于研究这些神经网络。需要使用NeuroLab库。
1. 创建一个感知器：
感知器是一个单独的神经元，负责执行所有的计算。这是一个简单的模型，但奠定了构建复杂神经网络的基础。该神经元将多个输入用不同的权重系数融合起来，并加上偏差值来计算输出。这是一个线性方程，它将输入值与感知器的输出关联起来。
import numpy as np
import neurolab as nl
import matplotlib.pyplot as plt
定义输入数据：
data=np.array([[0.3,0.2],[0.1,0.4],[0.4,0.6],[0.9,0.5]])
labels=np.array([[0],[0],[0],[1]])
画出点，以查看这些点的布局：
plt.figure()
plt.scatter(data[:,0],data[:,1])
plt.xlabel('X-axis')
plt.ylabel('Y-axis')
plt.title('Input data')
定义一个感知器perceptron，它有两个输入，该函数还需要限定输入数据的最大值和最小值：
perceptron=nl.net.newp([[0,1],[0,1]],1)
接下来训练感知器，epochs的数量指定了训练数据集需要完成的测试次数，show参数指定了显示训练过程的频率，lr参数指定了感知器的西、学习速度。学习速度是指学习算法在参数空间中搜索的步长，如果这个值太大，算法进行得会很快，但可能错失最优值；而如果这个值太小，算法可以命中最优值，但算法进行的会很慢，所以需要权衡，这里取0.01。
error=perceptron.train(data, labels, epochs=50, show=15, lr=0.01)
画出结果：
plt.figure()
plt.plot(error)
plt.xlabel('Number of epochs')
plt.ylabel('Training error')
plt.grid()
plt.title('Training error progress')
plt.show()
2. 创建一个单层神经网络：
单层神经网络由一个层次中的多个神经元组成。总体上看，单层神经网络将会有一个输入层、一个隐藏层和一个输出层。
import numpy as np
import matplotlib.pyplot as plt
import neurolab as nl
数据在文件中 需要加载：
input_file='data_single_layer.txt'
input_text=np.loadtxt(input_file)
data=input_text[:,0:2]
labels=input_text[:,2:]
画出输入数据：
plt.figure()
plt.scatter(data[:0],data[:,1])
plt.xlabel('X-axis')
plt.ylabel('Y-axis')
plt.title('Input data')
提取最小值和最大值：
x_min,x_max=data[:,0].min(), data[:,0].max()
y_min,y_max=data[:,1].min(), data[:,1].max()
定义一个单层神经网络，该网络的隐藏层包含两个神经元：
single_layer_net=nl.net.newp([[x_min, x_max], [y_min, y_max]], 2)
通过50次迭代训练神经网络：
error=single_layer_net.train(data, labels, epochs=50, show=20, lr=0.01
画出结果：
plt.figure()
plt.plot(error)
plt.xlabel('Number of epochs')
plt.ylabel('Training error')
plt.title('Training error progress')
plt.grid()
plt.show()
用新的测试数据来测试神经网络：
print(single_layer_net.sim([[0.3,4.5]]))
print(single_layer_net.sim([[4.5,0.5]]))
print(single_layer_net.sim([[4.3,8]]))
3. 创建一个深度神经网络：
创建的深度神经网络有一个输入层、多个隐藏层和一个输出层组成。
import neurolab as nl
import numpy as np
import matplotlib.pyplot as plt
定义参数，用于生成训练数据：
min_value=-12
max_value=12
num_datapoints=90
训练数据将由定义的一个函数组成，该函数将转换值，神经网络可以根据提供的输入数据和输出数据学习：
x=np.linspace(min_value, max_value, num_datapoints)
y=2*np.square(x)+7
y/=np.linalg.norm(y)
数组变形：
data=x.reshape(num_datapoint, 1)
labels=y.reshape(num_datapoints, 1)
画出输入数据：
plt.figure()
plt.scatter(data,labels)
plt.xlabel('X-axis')
plt.ylabel('Y-axis')
plt.title('Input data')
定义一个深度神经网络，该神经网络包含两个隐藏层，每个隐藏层包含10个神经元：
multilayer_net=nl.net.newff([[,om_value, max_value]],[10,10,1])
设置训练算法为梯度下降法：
multilayer_net.trainf=nl.train.train_gd
训练网络：
error=multilayer_net.train(data, labels, epochs=800, show=100, goal=0.01)
用训练数据运行该网络，查看其性能表现：
predicted_output=multilayer_net.sim(data)
画出训练误差结果：
plt.figure()
plt.plot(error)
plt.xlabel('Number of epochs')
plt.ylabel('Error')
plt.title('Training error progress')
创建一组新的输入数据，并运行神经网络，查看其性能表现：
x2=np.linspace(min_value, max_value, num_datapoints*2)
y2=multilayer_net.sim(x2.reshape(x2.size,1)).reshape(x2.size)
y3=predicted_output.reshape(num_datapoints)
画出输出结果：
plt.figure()
plt.plot(x2, y2, '-', x, y, '-', x, y3, 'p')
plt.title('Ground truth vs predicted output')
plt.show()
4. 创建一个向量量化器：
可以使用神经网络来做向量量化。向量量化是N维空间的舍入操作，广泛用于各个领域，如计算机视觉、自然语言处理和机器学习。
import numpy as np
import matplotlib.pyplot as plt
import neurolab as nl
数据在文件中 需要加载：
input_file='data_vq.txt'
input_text=np.loadtxt(input_file)
data=input_text[:,0:2]
labels=input_text[:,2:]
定义一个两层的学习向量量化LVQ（Learning Vector Quantization）神经网络，函数中的最后一个参数的数组指定了每个输出的加权百分比，各加权百分比之和应为1。
net=nl.net.newlvq(nl.tool.minmax(data), 10, [0.25, 0.25, 0.25, 0.25])
训练LVQ神经网络：
error=net.train(data, labels, epochs=100, goal=-1)
创建一个用于测试及可视化的网格点值：
xx, yy=np.meshgrid(np.arange(0, 8, 0.2), np.arange(0, 8, 0.2))
xx.shape=xx.size, 1
yy.shape=yy.size,1
input_grid=np.concatenate((xx,yy), axis=1)
用这些点值评价该网络：
output_grid=net.sim(input_grid)
在数据中定义4个类：
class1=data[labels[:,0]==1]
class2=data[labels[:,1]==1]
class3=data[labels[:,2]==1]
class4=data[labels[:,3]==1]
为每个类定义网络：
grid1=input_grid[output_grid[:,0]==1]
grid2=input_grid[output_grid[:,1]==1]
grid3=input_grid[output_grid[:,2]==1]
grid4=input_grid[output_grid[:,3]==1]
画出输出结果：
plt.plot(class1[:,0], class1[:,1], 'ko', class2[:,0], class2[:,1], 'ko', class3[:,0], class3[:,1], 'ko', class4[:,0], class4[:,1], 'ko')
plt.plot(grid1[:,0], grid1[:,1], 'b.', grid2[:,0], grid2[:,1], 'gx', grid3[:,0], grid3[:,1], 'cs', grid4[:,0], grid4[:,1], 'ro')
plt.axis([0,8,0,8])
plt.scatter(data[:0],data[:,1])
plt.xlabel('X-axis')
plt.ylabel('Y-axis')
plt.title('Vector quantization using neural networks')
plt.show()
5. 为序列数据分析创建一个递归神经网络：
递归神经网络能较好地分析序列和时间序列数据，数据的时序相关性非常关键。
import numpy as np
import matplotlib.pyplot as plt
import neurolab as nl
定义一个函数，该函数利用输入参数创建一个波形：
def create_waveform(num_points):
data1=1*np.cos(np.arange(0, num_points))
data2=2*np.cos(np.arange(0, num_points))
data3=3*np.cos(np.arange(0, num_points))
data4=4*np.cos(np.arange(0, num_points))
# 为每个区间创建不同的振幅，以此来创建一个随机波形
amp1=np.ones(num_points)
amp2=4+np.zeros(num_points)
amp3=2*np.ones(num_points)
amp4=0.5+np.zeros(num_points)
# 将数组合并生成输出数组，其中数据对应输入，而振幅对应相应的标签
data=np.array([data1, data2, data3, data4]).reshape(num_points*4, 1)
amplitude=np.array([amp1, amp2, amp3, amp4]).reshape(num_points*4, 1)
return data,amplitude
定义一个函数，用于画出将数据传入训练的神经网络后的输出：
def draw_output(net, num_points_test):
data_test, amplitude_test=create_waveform(num_points_test)
output_test=net.sim(data_test)
plt.plot(amplitude_test.reshape(num_points_test*4))
plt.plot(output_test.reshape(num_points_test*4))
定义mian函数，并生成示例数据：
if __name__=='__main__':
num_points=30
data, amplitude=create_waveform(num_points)
# 创建一个两层的递归神经网络
net=nl.net.newelm([[-2,2]], [10,1], [nl.trans.TanSig(), nl.trans.PureLin()])
# 设定每层的初始化函数
net.layers[0].initf=nl.init.InitRand([-0.1, 0.1], 'wb')
net.layers[1].initf=nl.init.InitRand([-0.1, 0.1], 'wb')
net.init()
# 训练递归神经网络
error=net.train(data, amplitude, epochs=1000, show=100, goal=0.01)
# 为训练数据计算来自网络的输出
output=net.sim(data)
# 画出训练误差
plt.subplot(211)
plt.plot(error)
plt.xlabel('Number of epochs')
plt.ylabels('Error (MSE)')
# 画出结果
plt.subplot(212)
plt.plot(amplitude.reshape(num_points*4))
plt.plot(output.reshape(num_points*4))
plt.legend(['Ground truth', 'Predicted output'])
# 创建一个随机长度的波形，查看该神经网络能否预测
plt.figure()
plt.subplot(211)
draw_output(net, 74)
plt.xlim([0, 300])
# 创建另一个长度更短的波形，查看该神经网络能否预测：
plt.subplot(212)
draw_output(net, 54)
plt.xlim([0, 300])
plt.show()
6. 在光学字符识别数据库中将字符可视化：
光学字符识别是指识别图像中的手写字符的过程。
import os
import sys
import cv2
import numpy as np
定义输入文件名：
input_fiel='letter.data'
定义可视化参数：
scaling_factor=10
start_index=6
end_index=-1
h,w=16,8
循环迭代文件直至用户按下Esc键，用Tab分隔符将行分隔成字符：
with open(input_file, 'r') as f:
for line in f.readlines():
data=np.array([255*float(x) for x in line.split('\t') [start_index:end_index]])
# 将数组重新调整为所需的形状，调整大小并将其展示
img=np.reshape(data, (h,w))
img_scaled=cv2.resize(img, None, fx=scaling_factor, fy=scaling_factor)
cv2.inshow('Image', img_scaled)
# 如果用户按下Esc键，则终止循环
c=cv2.waitKey()
if c==27:
break
7. 用神经网络创建一个光学字符识别器：
import numpy as np
import neurolab as nl
定义输入文件：
input_file='letter.data'
在用神经网络处理大量数据时，往往需要花费很多时间来做训练，这里只使用20个数据点：
num_datapoints=20
观察数据，可以看到在前20行有7个不同的字符。将其定义如下：
orig_labels='omandig' # 不同的字符
num_output=len(orig_labels) #不同字符的数量
用数据集的90%做训练，剩下的10%做测试，定义训练和测试参数为：
num_train=int(0.9*num_datapoints)
num_test=num_points-num_train
数据文件中每行的起始索引值和终止索引值设置为：
start_index=6
end_index=-1
生成数据集：
data=[]
labels=[]
with open(input_file, 'r') as f:
for line f.readlines():
list_vals=line.split('\t') # 按Tab键分割
# 如果字符不在标签列表中，跳过
if list_vals[1] not in orig_labels:
continue
# 提取标签，并将其添加到主列表的后面
label=np.zeros((num_output, 1))
label[orig_labels.index(list_vals[1])]=1
labels.append(label)
# 提取字符，并将其添加到主列表的后面
cur_char=np.array([float(x) for x in list_vals[start_index:end_index]])
data.append(cur_char)
# 当有足够多数据时跳出循环
if len(data)>=num_datapoints:
break
将以上数据转换成NumPy数组：
data=np.asfarray(data)
labels=np.array(labels).reshape(num_datapoints, num_output)
提取数据的维度信息：
num_dims=len(data[0])
用10000次迭代来训练神经网络：
net=nl.net.newff([[0,1] for _ in range(len(data[0]))], [128,16,num_output])
net.trainf=nl.train.train_gd
error=net.train(data[:num_train,:], labels[:num_train,:], epochs=10000, show=100, goal=0.01)
为测试输入数据预测输出结构：
predicted_output=net.sim(data[num_train:,:])
for i in range(num_test):
print('\nOriginal:', orig_labels[np.argmax(labels[i])]
print('Predicted:', orig_labels[np.argmax(predicted_output[i])]
十三、可视化数据：
数据的视觉可以帮助选择正确的算法，数据可视化的主要目标之一就是用图和表清晰地表达出数据，以便更准确、更有效地交流信息。用不同的图来展示各个变量之间的模式或关系。
1. 画3D散点图：
import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
生成一个空白图像：
fig=plt.figure()
ax=fig.add_subplot(111, projection='3d')
定义生成的值的个数：
n=250
生成一个lambda函数来生成给定范围的值：
f=lambda minval,maxval,n:minval+(maxval-minval)*np.random.rand(n)
用这个函数生成X、Y和Z值：
x_vals=f(15, 41, n)
y_vals=f(-10, 70, n)
x_vals=f(-52, -37, n)
画出这些值：
ax.scatter(x_vals, y_vals, z_vals, c='k', marker='o')
ax.set_xlabel('X axis')
ax.set_ylabel('Y axis')
ax.set_zlabel('Z axis')
plt.show()
2. 画气泡图：
import numpy as np
import matplotlib.pyplot as plt
定义生成的值的个数：
n=40
生成随机的x值和y值：
x=np.random.rand(num_vals)
y=np.random.rand(num_vals)
在气泡图中定义每个点的面积值：
max_radius=25
area=np.pi*(max_radius*np.random.rand(num_vals))**2
定义颜色：
colors=np.random.rand(num_vals)
画出这些值：
plt.scatter(x, y, s=area, c=colors, alpha=1.0)
plt.show()
3. 画动态气泡图：
import numpy as np
import matplotlib.pyplot as plt
from matplotlib.animation import FuncAnimation
定义一个tracker函数，该函数将动态更新气泡图：
def tracker(cur_num):
cur_index=cur_num%num_points
定义颜色：
datapoints['color'][:,3]=1.0
更新圆圈的大小：
datapoints['size']+=datapoints['growth']
更新集合中最老的数据点的位置：
datapoints['position'][cur_index]=np.random.uniform(0,1,2)
datapoints['size'][cur_index]=7
datapoints['color'][cur_index]=(0,0,0,1)
datapoints['growth'][cur_index]=np.random.uniform(40,150)
更新散点图的参数：
scatter_plot.set_edgecolors(datapoints['color'])
scatter_plot.set_sizes(datapoints['size'])
scatter_plot.set_offsets(datapoints['position'])
定义一个主函数并生成一个空白图像：
if __name__=='__main__':
fig=plt.figure(figsize=(9,7), facecolor=(0, 0.9, 0.9))
ax=fig.add_axes([0, 0, 1, 1], frameon=False)
ax.set_xlim(0,1), ax.set_xticks([])
ax.set_ylim(0,1), ax.set_yticks([])
定义在任意时间点上的点的个数：
num_points=20
用随机值定义这些数据点：
datapoints=np.zeros(num_points, dtype=[('position', float, 2), ('size', float, 1), ('growth', float, 1), ('color', float, 4)])
datapoints['position']=np.random.uniform(0, 1, (num_points, 2))
datapoints['growth']=np.random.uniform(40, 150, num_points)
创建一个散点图，该散点图的每一帧都会更新：
scatter_plot=ax.scatter(datapoints['position'][:,0], datapoints['position'][:,1], s=datapoints['size'], lw=0.7, edgecolors=datapoints['color'], facecolors='none')
用tracker函数启动动态模拟：
animation=FuncAnimation(fig, tracker, interval=10)
plt.show()
4. 画饼图：
import numpy as np
import matplotlib.pyplot as plt
定义各标签和相应的值：
data={'Apple':26, 'Mango':17, 'Pineapple':21, 'Banana':29, 'Strawberry':11}
定义可视化颜色：
colors=['orange', 'lightgreen', 'lightblue', 'gold', 'cyan']
定义一个变量，以突出饼图的一部分，将其与其他部分分离开。如果不想突出任何部分，将所有的值设置为0：
explode=(0,0.2,0,0,0)
画饼图：
plt.pie(list(data,values(), explode, labels=data.keys(), colors=colors, autoopct='%1.1f%%', shadow=False, starangle=90)
设置饼图的宽高比，equal表示为圆形：
plt.axis('equal')
plt.show()
5. 画日期格式的时间序列数据：
import numpy as np
import matplotlib.pyplot as plt
from matplotlib.nlab import csv2rec
import matplotlib.cbook as cbook
firm matplotlib.ticker import Formatter
定义一个用于将日期格式化的类。init函数设置类变量：
class DataFormatter(Formatter):
def __init__(self, dates , date_format='%Y-%m-%d'):
self.dates=dates
self.date_format=date_format
# 提取给定时间的值，并用如下格式将其返回：
def __call__(self, t, position=0):
index=int(round(t))
if index>=len(self.dates) or index<0:
return ''
return self.dates[index].strftime(self.date_format)
定义main函数，需要用到matplotlib中苹果公司的股票报价CSV文件：
if __name__='__main__':
input_file=cbook.get_sample_data('aapl.csv', asfileobj=False)
# 加载CSV文件中数据到NumPy数组
data=csv2rec(input_file)
# 提取这些值的子集
data=data[-70:]
# 创建一个格式化对象，并将其用日期数据初始化
formatter=DateFormatter(data.date)
定义X轴和Y轴
x_vals=numpy.arange(len(data))
y_vals=data.close # y轴表示收盘价
# 画出数据：
fig,ax=plt.subplots()
ax.xaxis.set_major_formatter(formatter)
ax.plot(x_vals, y_vals, 'o-')
fig.autofmt_xdate()
plt.show()
6. 画直方图：
如果需要对比两组数据，可以使用直方图做对比。
import numpy as np
import matplotlib.pyplot as plt
本例对比苹果和橘子的产量，数据为：
apples=[30, 25, 22, 36, 21, 29]
oranges=[24, 33, 19, 27, 35, 20]
num_groups=len(apples)
创建一个图像并定义其参数：
fig,ax=plt.subplots()
定义x轴：
indices=np.arange(num_groups)
直方图的宽度和透明度：
bar_width=0.4
opacity=0.6
画直方图：
hist_apples=plt.bar(indices, apples, bar_width, alpha-opacity, color='g', label='Apples')
hist_oranges=plt.bar(indices+bar_width, oranges, bar_width, alpha-opacity, color='b', label='Oranges')
设置直方图的参数：
plt.xlabel('Month')
plt.ylabel('Production quantity')
plt.title('Comparing apples and oranges')
plt.xticks(indices+bar_width, ('Jan', 'Feb', 'Mar', 'Apr', 'May', 'Jun'))
plt.ylim([0,45])
plt.legend()
plt.tight_layout()
plt.show()
7. 可视化热力图：
当两组数据的各点具有一定的相关性时，可以用这种图形表现方式。矩阵包含的单个值都被表示为图中的颜色值。
import numpy as np
import matplotlib.pyplot as plt
定义两组数据：
group1=['France', 'Italy', 'Spain', 'Portugal', 'Germany']
group2=['Japan', 'China', 'Brazil', 'Russia', 'Australia']
生成一个随机二维矩阵：
data=np.random.rand(5,5)
创建一个图像：
fig,ax=plt.subplots()
创建一个热力图：
heatmap=ax.pcolor(data, cmap=plt.cm.gray)
画出这些值：
ax.set_xticks(np.arange(data.shape[0])+0.5, minor=False)
ax.set_yticks(np.arange(data.shape[1])+0.5, minor=False)
ax.invert_yaxis()
ax.xaxis.tick_top()
ax.set_xticklabels(group2, minor=False0
ax.set_yticklabels(group1, minor=False0
plt.show()
8. 动态信号的可视化模拟：
如果需要将实时信号可视化，最好的方式是看到波形的产生。这里是对实时动态变化的信号进行模拟并可视化。
import numpy as np
import matplotlib.pyplot as plt
import matplotlib.animation as animation
创建一个函数，用于生成阻尼正弦信号：
def generate_data(length=2500, t=0, step_size=0.05):
for count in range(length):
t+=step_size
signal=np.sin(2*np.pi*t)
damper=np.exp(-t/8.0)
yield t, signal*damper
定义一个initializer函数，用于将图像中的参数初始化：
def initializer():
peak_val=1.0
buffer_val=0.1
# 设置这些参数
ax.set_ylim(-peak_val*(1+buffer_val), peak_val*(1+buffer_val))
ax.set_xlim(0,10)
del x_vals[:]
del y_vals[:]
line.set_data(x_vals, y_vals)
return line
定义一个函数画出数值：
def draw(data):
t,signal=data
x_vals.append(t)
y_vals.append(signal)
x_min,x_max=ax.get_xlim()
# 如果这些值超出当前X轴最大值范围，就更新X轴最大值并扩展图像：
if t>=x_max:
ax.set_xlim(x_min, 2*x_max)
ax.figure.canvas.draw()
line.set_data(x_vals, y_vals)
return line
定义主函数：
if __name__='__main__':
fig.ax=plt.subplots()
ax.grid()
# 提取线
line,=ax.plot([],[],lw=1.5)
# 创建变量，并用空列表对其初始化
x_vals, y_vals=[],[]
# 用动画器对象定义并启动动画
animator=animation.FuncAnimation(fig, draw, generate_data, blit=False, interval=10, repeat=False, init_func=initializer)
plt.show()
1069
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
