第一个极小的机器学习的应用

Posted 2020-08-15 曹孟德

  现在给出一个Web统计信息，他们存储着每小时的访问次数。每一行包含连续的小时和信息，以及该小时Web的访问次数。现在要解决的问题是，估计在何时访问量达到基础设施的极限。极限数据是每小时100000次访问。

1.读取数据：

# 获取数据
filepath = r\'C:\\Users\\TD\\Desktop\\data\\Machine Learning\\1400OS_01_Codes\\data\\web_traffic.tsv\'
data = sp.genfromtxt(filepath,delimiter = \'\\t\')
x = data[:,0]
y = data[:,1]

其中，x表示小时，y表示访问量。

2.预处理和清洗数据：

print sp.sum(sp.isnan(y))

结果显示含有8个控值，为了方便，在此处理缺失值办法是直接剔除。

x = x[~sp.isnan(y)]
y = y[~sp.isnan(y)]

接下来，画出散点图，观察数据的规律：

# 可视化，观察数据规律
plt.scatter(x,y)
plt.title(\'Web traffic over the last month\')
plt.xlabel(\'Time\')
plt.ylabel(\'Hits/hours\')
plt.xticks([w*24*7 for w in range(5)],
		   [\'week {}\'.format(i) for i in range(5)])
plt.autoscale(tight = True)
plt.grid()
plt.show()

3 选择正确的模型和学习算法：

回答原始问题需要明确以下几点：

1）找到噪声数据后真正的模型

2）使用这个模型预测未来，一遍解决我们的问题

1.首先需要明白模型与实际数据区别。模型可以理解为对复杂现实世界简化的理论近似。它总会包含一些劣质的类容，这个就叫做近似误差。我们用真实数据与模型预测的数据之间的平方距离来计算这个误差，对于一个训练好的模型f，按照下面函数来计算误差：

def error(f,x,y):
	return sp.sum((f(x)-y)**2)

 2.简单的线性模型

现在用一个线性模型来拟合上面数据，看看可以得到什么。

plt.scatter(x,y)
plt.title(\'Web traffic over the last month\')
plt.xlabel(\'Time\')
plt.ylabel(\'Hits/hours\')
plt.xticks([w*24*7 for w in range(5)],
		   [\'week {}\'.format(i) for i in range(5)])
plt.autoscale(tight = True)
plt.grid()

# 开始构建模型,使用1阶多项式拟合
p1 = sp.polyfit(x,y,1)
f1 = sp.poly1d(p1)  # 将拟合系数传入ployld函数创建模型函数f1
fx = sp.linspace(0,x[-1],1000)
plt.plot(fx,f1(fx),linewidth = 3)
plt.legend(["d = {}".format(f1.order)], loc = "upper left")
plt.show()

　

上图显示了第一个训练的模型，发现前四个星期好像没有偏差很多，可以清楚的看到直线模型的假设是有问题的。

3 接下来用3阶，10阶，50阶多项式来拟合：

colors = [\'g\', \'k\', \'b\', \'m\', \'r\']

def error(f,x,y):
	return sp.sum((f(x)-y)**2)

def plot_models(x, y, models, fname, mx = None):
    plt.clf()
    plt.scatter(x, y, s=10)
    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)])

    if models:
        if mx is None:
            mx = sp.linspace(0, x[-1], 1000)
        for model, color in zip(models, colors):
            # print "Model:",model
            # print "Coeffs:",model.coeffs
            plt.plot(mx, model(mx), c = color,linewidth = 1.5)

        plt.legend(["d = {}".format(m.order) for m in models], loc="upper left")
    plt.autoscale(tight=True)
    plt.grid(True, linestyle=\'-\', color=\'0.75\')
    plt.savefig(fname)

# create and plot models
os.chdir(r\'C:\\Users\\TD\\Desktop\\data\\Machine Learning\\1400OS_01_Codes\\data\')
f1 = sp.poly1d(sp.polyfit(x,y,1))
f3 = sp.poly1d(sp.polyfit(x, y, 3))
f10 = sp.poly1d(sp.polyfit(x, y, 10))
f50 = sp.poly1d(sp.polyfit(x, y, 50))
plot_models(
    x, y, [f1,f3, f10, f50],"2.png")

# error
indices = [1,3,10,50]
for index,model in zip(indices,[f1,f3,f10,f50]):
	print \'Error d= {} : {}\'.format(index,error(model,x,y))

可以看出多项式越复杂，数据逼近越好。他们误差如下：

Error d= 1 : 317389767.34
Error d= 3 : 139350144.032
Error d= 10 : 121942326.364
Error d= 50 : 109504607.366

看看10阶和100阶的多项式，我们发现了巨大的震荡。似乎这样拟合的太过了，他不断捕捉到背后数据的生成，还把噪声数据也考虑进去了。这样叫做过拟合。然而1阶的显然太简单了，不能反映数据的规律，这种叫做欠拟合。不管是欠拟合还是过拟合，都不适合进行预测。

4 已退为进，另眼看数据

观察数据，似乎第三周和第四周之间有一个拐点。这可以让我们将3.5周作为分界点，把数据分为两份，并训练出两条直线。

plt.scatter(x, y, s=10)
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)])
inflection = int(3.5*7*24)
xa = x[:inflection]
ya = y[:inflection]
xb = x[inflection:]
yb = y[inflection:]
fa = sp.poly1d(sp.polyfit(xa,ya,1))
fb = sp.poly1d(sp.polyfit(xb,yb,1))
plt.scatter(x,y)
fax = sp.linspace(0,x[-1],1000)
fbx = sp.linspace(x[inflection]/1.1,x[-1]*1.1,1000)
plt.plot(fax,fa(fax),c = \'g\',linewidth = 2.5)
plt.plot(fbx,fb(fbx),c = \'r\',linewidth = 2.5)
plt.show()

　

很明显，这两条直线组合起来似乎比之前任何模型都可以更好的拟合数据，虽然组合后的误差高于高阶多项式的误差。为什么仅仅在最后一周上更相信线性模型呢？这是因为我们认为他更好的符合未来数据。10阶和100阶多项式在此没有光明的未来，他们只是非常努力的对给定的数据进行拟合，但是他们却无法推广到将来的数据上，这就是过拟合，另外低阶模型也不能恰好的模拟数据，叫做欠拟合。

5 训练与测试

如果有些外来数据用于模型评估，那么仅从近似误差结果就可以判断出我们的选择的模型是好还是坏。尽管我们找不到未来数据，但是可以从现有的数据中拿出一部分，来判断我们的结果是好还是坏了。利用拐点后的数据进行训练，得到的二阶模型的误差最小，这个模型很适中，既不欠拟合也不过拟合。

6 回答最初的问题

得到了训练的模型。只需要带入数值就可以计算得到我们所求结果。

实验代码：

#!/usr/bin/env python
# -*- coding: utf-8 -*-
# __author__ : \'小糖果\'

import scipy as sp
import matplotlib.pyplot as plt
from scipy.optimize import fsolve
import os

# 获取数据
filepath = r\'C:\\Users\\TD\\Desktop\\data\\Machine Learning\\1400OS_01_Codes\\data\'
data = sp.genfromtxt(os.path.join(filepath,\'web_traffic.tsv\'),delimiter = \'\\t\')
x = data[:,0]
y = data[:,1]

# 缺失数据处理，用相邻数据平均数代替
print sp.sum(sp.isnan(y))
x = x[~sp.isnan(y)]
y = y[~sp.isnan(y)]

colors = [\'g\', \'k\', \'b\', \'m\', \'r\']

def error(f,x,y):
	return sp.sum((f(x)-y)**2)

def plot_model(x, y,models = None,fname = None,mx = None):
	plt.clf()
	plt.scatter(x,y)
	plt.title("Web traffic over the last month")
	plt.xlabel("Time")
	plt.ylabel("Hits/hour")
	plt.xticks([w*24*7 for w in range(10)],
			   [\'week {}\'.format(i) for i in range(10)])
	if models:
		if mx is None:
			mx = sp.linspace(0,x[-1],1000)
		for (model,color) in zip(models,colors):
			plt.plot(mx,model(mx),c = color,linewidth = 2)
		plt.legend([\'d = {}\'.format(m.order) for m in models],loc = \'upper left\')
	plt.autoscale(tight = True)
	plt.grid(True)
	if fname:
		plt.savefig(fname)
	else:
		plt.show()

# 查看初始数据
plot_model(x,y,fname = os.path.join(filepath,\'1.jpg\'))

#分别用1，2,10,50阶多项式拟合
f1 = sp.poly1d(sp.polyfit(x,y,1))
f2 = sp.poly1d(sp.polyfit(x,y,2))
f10 = sp.poly1d(sp.polyfit(x,y,10))
f50 = sp.poly1d(sp.polyfit(x,y,50))
plot_model(x,y,models = [f1,f2,f10,f50],fname = os.path.join(filepath,\'2.jpg\'))

# 线性分段拟合
plt.clf()
inflection = int(3.5*7*24)
xa = x[:inflection]
ya = y[:inflection]
xb = x[inflection:]
yb = y[inflection:]
fa = sp.poly1d(sp.polyfit(xa,ya,1))
fb = sp.poly1d(sp.polyfit(xb,yb,1))
plt.scatter(x,y)
fax = sp.linspace(0,x[-1],1000)
fbx = sp.linspace(x[inflection]/1.1,x[-1]*1.1,1000)
plt.plot(fax,fa(fax),c = \'g\',linewidth = 2.5)
plt.plot(fbx,fb(fbx),c = \'r\',linewidth = 2.5)
plt.grid(True)
plt.savefig(os.path.join(filepath,\'3.jpg\'))

# 只使用后面部分数据
f1 = sp.poly1d(sp.polyfit(xb,yb,1))
f2 = sp.poly1d(sp.polyfit(xb,yb,2))
f10 = sp.poly1d(sp.polyfit(xb,yb,10))
f50 = sp.poly1d(sp.polyfit(xb,yb,50))
plot_model(xb,yb,models = [f1,f2,f10,f50],
		   mx = sp.linspace(xb[0],xb[-1],100),
		   fname = os.path.join(filepath,\'4,jpg\'))

# 求问题的解,使用二次多项式模型
ans = fsolve(f2 - 100000,800)/7/24
print ans