# 获取数据
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]

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

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()

def error(f,x,y):return sp.sum((f(x)-y)**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()

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.coeffsplt.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))

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()

#!/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