building machine learning system with Python 学习笔记--从零开始机器学习（3）第一个应用

这个小应用是根据已有的网站访问量来预测什么时候到达现有设施的极限，估计是每小时100000个请求。

这种问题在初高中求函数极值时经常遇到，只是现在函数形式是未知的，只有一定量的离散数据。机器学习就派上用场了，通过猜测函数形式，估算损失函数，看哪种函数估算时错误最少。

查看已有数据：该数据文件包含在building machine learning system with Python 学习笔记--从零开始机器学习（0）中文学习资料中

读取文件数据，打印出了前十行

>>> import scipy as sp

有743个二维数据点

1）预处理和清洗数据

python -i 可以在交互环境中使用文件，x是第一列，y是第二列，sp.sum统计y中非法数据有8个

x = x[~sp.isnan(y)]对非法数据的下标取反，过滤掉非法数据。并使用matplotlib库将散点图画出来（在http://matplotlib.org/users/pyplot_tutorial.html上看更多画图教程）：

2）选择正确的模型和学习算法

a)找到有噪数据背后真正的模型

b)使用该模型预测

构建第一个模型：

准备工作--定义误差函数：

def error(f, x, y):

     return sp.sum((f(x)-y)**2)

选择简单直线作为模型：

polyfit()函数会把拟合的模型函数所使用的参数返回，即fp1，full设置为True可以获得更多逼近过程的背景信息

最优的近似直线如下：

f(x) = 2.59619213 * x + 989.02487106,

残差是3.17389767*(10**8)

用poly1d()根据这些参数创建一个模型函数

到这一步的代码如下，fx=sp.linspace生成x的值，linspace将从0到x最后一个元素的值分为1000等分

效果如下：d=1表示f1是一维函数

没有比较不知道直线的残差3.17389767*(10**8)是好还是坏，将直线设为基准，把模型设置成二维的

为了便于对比，修改后的代码

展示：红色虚线是二维模型

误差179983507.878几乎是直线模型误差的一半

近似多项式：f(x) = 0.0105322215 * x**2 - 5.26545650 * x + 1974.76082

这样看来，复杂性是不是越大越好?不妨将阶数设为3、 10、 100

展示：

polyfit维度设置为100时，展示的却是53，原因：http://stackoverflow.com/questions/20838970/scipy-polyfitx-y-100-would-be-100th-order-polynome-but-matplotlib-pyplot-le

In fact, most significant component seems to have degree 2. So it's normal, that best approximating your data polinomial of degree not greater than 100, in fact has degree 53.All higher monomials are degenerate.

有人的回答是53就是最优拟合了，再高阶实际上是在退化。anyway, I can not figure it out myself。有懂的同学还请指教。

误差结果：10阶和100阶的多项式不但捕捉到了背后的数据生成过程，还把噪声也包含进去了，这叫过度拟合（overfitting）

现有的选择：

选择其中一个拟合出的多项式模型

换成另外一类更复杂的模型

从不同的角度思考数据，然后重新开始

在5个拟合模型中1阶太简单，10阶和100阶过度拟合，2阶和3阶似乎比较匹配，但是进行预测时效果并不好。

重新审视数据：

第三周和第四周的数据之间有一个拐点，可以尝试以3.5周为分界点，训练出两条直线来

import scipy as spdata = sp.genfromtxt("web_traffic.tsv", delimiter="\t")x = data[:,0]y = data[:,1]x = x[~sp.isnan(y)]y = y[~sp.isnan(y)]inflection = int(3.5*7*24)xa = x[:inflection]ya = y[:inflection]xb = x[inflection:]yb = y[inflection:]def error(f, x, y):    return sp.sum((f(x)-y)**2)fa = sp.poly1d(sp.polyfit(xa, ya, 1))fb = sp.poly1d(sp.polyfit(xb, yb, 1))fa_error = error(fa, xa, ya)fb_error = error(fb, xb, yb)print("Error inflection=%f" % (fa_error + fb_error))import matplotlib.pyplot as pltfig = plt.figure(figsize=(5,5))ax = fig.add_subplot(111)plt.scatter(x, y)plt.title("Web traffic over the last month")plt.xlabel("Time")plt.ylabel("Hits/hour")plt.xticks([w*7*24 for w in range(10)], ['week %i' %w for w in range(10)])plt.autoscale(tight=True)fxa = sp.linspace(0, xa[-1], 1000)fxb = sp.linspace(xa[-1], xb[-1], 1000)ax.plot(fxa, fa(fxa), 'g--', linewidth=4,label='d=%i'% fa.order)ax.plot(fxb, fb(fxb), 'r--', linewidth=4, label='d=%i'% fb.order)plt.grid()plt.show()

展示：

打印残差：中文版有一处错误，print("Error inflection=%f" %(fa + fb_error)) fa是poly1d对象，和float数值相加会报错

可以看到两条线组合起来似乎比之前的模型都能更好地拟合数据。但组合之后的误差仍然高于高阶多项式的误差。

10阶和100阶的模型非常努力地对给定数据正确建模，但无法推广到将来的数据，叫做过度拟合（overfitting）。低阶模型似乎也不能恰当地拟合数据，叫欠拟合（underfitting）

只用拐点之后的数据来训练模型

import scipy as spdata = sp.genfromtxt("web_traffic.tsv", delimiter="\t")#print(data[:10])#print(data.shape)x = data[:,0]y = data[:,1]x = x[~sp.isnan(y)]y = y[~sp.isnan(y)]inflection = int(3.5*7*24)xb = x[inflection:]yb = y[inflection:]def error(f, x, y):    return sp.sum((f(x)-y)**2)fp1, res, rank, sv, rcond = sp.polyfit(xb, yb, 1, full=True)f1 = sp.poly1d(fp1)f2p = sp.polyfit(xb, yb, 2)f2 = sp.poly1d(f2p)f3p = sp.polyfit(xb, yb, 3)f3 = sp.poly1d(f3p)f4p = sp.polyfit(xb, yb, 10)f4 = sp.poly1d(f4p)f5p = sp.polyfit(xb, yb, 100)f5 = sp.poly1d(f5p)import matplotlib.pyplot as pltfig = plt.figure(figsize=(5,5))ax = fig.add_subplot(111)plt.scatter(xb, yb)plt.title("Web traffic over the last month")plt.xlabel("Time")plt.ylabel("Hits/hour")plt.xticks([w*7*24 for w in range(10)], ['week %i' %w for w in range(10)])plt.autoscale(tight=True)fx = sp.linspace(inflection, xb[-1], 100)ax.plot(fx, f1(fx), linewidth=4,label='d=%i'% f1.order)ax.plot(fx, f2(fx), 'r--', linewidth=4, label='d=%i'% f2.order)ax.plot(fx, f3(fx), 'g--', linewidth=4, label='d=%i'% f3.order)ax.plot(fx, f4(fx), 'y--', linewidth=4, label='d=%i'% f4.order)ax.plot(fx, f5(fx), 'b--', linewidth=4, label='d=%i'% f5.order)plt.legend(loc="upper left")plt.grid()plt.show()

展示