机器学习–线性回归
所谓回归分析(regression analysis),是确定两种或两种以上变量间相互依赖的定量关系的一种统计分析方法,也就是研究当自变量x变化时,因变量y以何种形式在变化。在机器学习领域,回归应用于被预测对象具有连续值特征的情况(如客流量、降雨量、销售量等)。
在机器学习的线性回归分析中,如果只包括一个自变量(特征x)和一个因变量(标签y),且两者的关系可用一条直线近似表示,这种回归分析就称为一元线性回归分析。如果回归分析中包括两个或两个以上的自变量,且因变量和自变量之间是线性关系,则称为多元线性回归分析。
举一个栗子
通过机器学习算法,根据过去记录下来的广告投放金额和商品销售额,来预测在未来的某个节点,一个特定的广告投放金额对应能实现的商品销售额。
(1)明确定义所要解决的问题----网店销售额的预测。
(2)在数据的收集和预处理环节,分5个小节完成数据的预处理工作,分别如下。
■将收集到的数据可视化,显示出来看一看。
■做特征工程,使数据更容易被机器处理。
■拆分数据集为训练集和测试集。
■做特征缩放,把数据值压缩到比较小的区间。
(3)选择机器学习模型的环节,其中有3个主要内容。
■确定机器学习的算法—这里也就是线性回归算法。
■确定线性回归算法的假设函数。
■确定线性回归算法的损失函数。
(4)通过梯度下降训练机器,确定模型内部参数的过程。
(5)进行超参数调试和性能优化。
数据的收集和预处理
收集网店销销售额数据
数据可视化
上传的数据要读取注意要填写好位置
数据相关分析
导入数据可视化所需要的库
相关分析(correlation analysis):正值表示正相关,负值表示负相关。数值越大,相关性越强。
如果a和b的相关性系数是1,则a和b总是相等的。
import matplotlib.pyplot as plt #Matplotlib – Python画图工具库
import seaborn as sns #Seaborn – 统计学数据可视化工具库
#对所有的标签和特征两两显示其相关性的热力图(heatmap)
sns.heatmap(df_ads.corr(), cmap="YlGnBu", annot = True)
plt.show() #plt代表英文plot,就是画图的意思
初步分析可得微信公众号里面做广告是最为合理的选择。
数据的散点图
#显示销量和各种广告投放量的散点图
sns.pairplot(df_ads,
x_vars=['wechat', 'weibo', 'others'],
y_vars='sales',
height=4, aspect=1, kind='scatter')
plt.show()
数据集清洗和规范化
在本案例的3个特征中,微信广告投放金额和商品销售额的相关性比较高。
为了简化模型,我们只留下微信广告投放金额数据。这样,就把多变量的回归分析简化为单变量的回归分析。
X = np.array(df_ads.wechat) #构建特征集,只含有微信广告一个特征
y = np.array(df_ads.sales) #构建标签集,销售金额
print ("张量X的阶:",X.ndim)
print ("张量X的形状:", X.shape)
print ("张量X的内容:", X)
目前X数组中只有一个特征,张量的阶为1,这个1D的特征张量不是机器学习算法能够接受的格式
对于回归问题的数值类型数据集,机器学习模型所读入的规范格式应该是2D张量,也就是矩阵,其形状为 (样本数,标签数)。
那么就现在的特征张量X而言,则是要把它的形状从(200,)变成(200,1),然后再进行机器学习。
因此需要用reshape方法给上面的张量变形:
X = X.reshape((len(X),1)) #通过reshape函数把向量转换为矩阵,len函数返回样本个数
y = y.reshape((len(y),1)) #通过reshape函数把向量转换为矩阵,len函数返回样本个数
print ("张量X的阶:",X.ndim)
print ("张量X的形状:", X.shape)
print ("张量X的内容:", X)
拆分为训练集和测试集
#将数据集进行80%(训练集)和20%(测试集)的分割
from sklearn.model_selection import train_test_split
X_train, X_test, y_train, y_test = train_test_split(X, y,
test_size=0.2, random_state=0)
Sklearn中的train_test_split
函数,是机器学习中拆分数据集的常用工具。
test_size=0.2,表示拆分出来的测试集占总样本量的20%。
如果用print语句输出拆分之后的新数据集(如X_train、X_test)的内容,会发现这个工具已经为数据集进行了乱序(重新随机排序)的工作,因为其中的shuffle参数默认值为True。
而其中的random_state
参数,则用于数据集拆分过程的随机化设
定。如果指定了一个整数,那么这个数叫作随机化种子,每次设定固定
的种子能够保证得到同样的训练集和测试集,否则进行随机分割。
数据归一化
归一化是按比例的线性缩放。数据分布不变,但是都落入一个小的特定区间,比如0~1或者-1~+1。
x ′ = x − min ( x ) max ( x ) − min ( x ) x^{\prime}=\frac{x-\min(x)}{\max(x)-\min(x)} x′=max(x)−min(x)x−min(x)
def scaler(train, test): #定义归一化函数,进行数据压缩
min = train.min(axis=0) #训练集最小值
max = train.max(axis=0) #训练集最大值
gap = max - min #最大值和最小值的差
train -= min #所有数据减最小值
train /= gap #所有数据除以大小值差
test -= min #把训练集最小值应用于测试集
test /= gap #把训练集大小值差应用于测试集
return train, test #返回压缩后的数据
X_train,X_test = scaler(X_train,X_test) #对特征归一化
y_train,y_test = scaler(y_train,y_test) #对标签也归一化
#用之前已经导入的matplotlib.pyplot中的plot方法显示散点图
plt.plot(X_train,y_train,'r.', label='Training data')
plt.xlabel('Wechat Ads') # x轴Label
plt.ylabel('Sales') # y轴Label
plt.legend() # 显示图例
plt.show() # 显示绘图结果
选择机器学习模型
图中的一元线性特征很明显,所以用一元线性函数就可以进行拟合。
如果模型的预测完全准确,则损失(loss)为0;如果不准确,就有损失。
根据loss创建的函数也叫损失函数。
在回归中,损失函数经常是:均方误差,平均绝对误差,平均偏差误差等
在分类中,损失函数经常是:交叉熵损失,多分类SVM损失,等
在这个工程中,我们用MES函数,即最小二乘法
L ( w , b ) = M S E = 1 2 N ∑ ( x , y ) ∈ D ( y − h ( x ) ) 2 L(w,b)=MSE=\frac{1}{2N}\sum_{(x,y)\in D}(y-h(x))^2 L(w,b)=MSE=2N1(x,y)∈D∑(y−h(x))2
def loss_function(X, y, weight, bias): # 手工定义一个MSE均方误差函数
y_hat = weight*X + bias # 这是假设函数,其中已经应用了Python的广播功能
loss = y_hat-y # 求出每一个y’和训练集中真实的y之间的差异
cost = np.sum(loss**2)/(2*len(X)) # 这是均方误差函数的代码实现
return cost # 返回当前模型的均方误差值
print ("当权重5,偏置3时,损失为:",
loss_function(X_train, y_train, weight=5, bias=3))
print ("当权重100,偏置1时,损失为:",
loss_function(X_train, y_train, weight=100, bias=1))
当权重5,偏置3时,损失为: 12.79639097078006
当权重100,偏置1时,损失为: 1577.9592615030556
通过梯度下降找最佳参数
梯度 = ∂ ∂ w L ( w ) = ∂ ∂ w 1 2 N ∑ ( x , y ) ∈ D ( y − h ( x ) ) 2 = 1 2 N ∑ ( x , y ) ∈ D ( y − ( w ∙ x ) ) ∙ x \text{梯度}=\frac{\partial}{\partial w}L(w)=\frac{\partial}{\partial w}\frac{1}{2N}\sum_{(x,y)\in D}(y-h(x))^2=\frac{1}{2N}\sum_{(x,y)\in D}(y-(w\bullet x))\bullet x 梯度=∂w∂L(w)=∂w∂2N1(x,y)∈D∑(y−h(x))2=2N1(x,y)∈D∑(y−(w∙x))∙x
即对权重求偏导
学习速率
w随梯度下降的更新如下
w = w − α ⋅ ∂ ∂ w L ( w ) w=w-\alpha\cdot\frac{\partial}{\partial w}L(w) w=w−α⋅∂w∂L(w)
在机器学习刚刚开始的时候,学习速率设置得大一些,当接近最佳权重时,减小学习速率。
def gradient_descent(X, y, w, b, lr, iter): # 定义一个实现梯度下降的函数
l_history = np.zeros(iter) # 初始化记录梯度下降过程中损失的数组
w_history = np.zeros(iter) # 初始化记录梯度下降过程中权重的数组
b_history = np.zeros(iter) # 初始化记录梯度下降过程中偏置的数组
for i in range(iter): # 进行梯度下降的迭代,就是下多少级台阶
y_hat = w*X + b # 这个是向量化运行实现的假设函数
loss = y_hat-y # 这是中间过程,求得的是假设函数预测的y和真正的y值间的差值
derivative_w = X.T.dot(loss)/len(X) # 对权重求导, len(X)是样本总数
derivative_b = sum(loss)*1/len(X) # 对偏置求导
w = w - lr*derivative_w # 结合下降速率alpha更新权重
b = b - lr*derivative_b # 结合下降速率alpha更新偏置
l_history[i] = loss_function(X, y, w,b) # 梯度下降过程中损失的历史
w_history[i] = w # 梯度下降过程中权重的历史
b_history[i] = b # 梯度下降过程中偏置的历史
return l_history, w_history, b_history # 返回梯度下降过程数据
实现一元线性回归模型并调试参数
# 首先确定参数的初始值
iterations = 100; # 迭代100次
alpha = 0.5; # 此处初始学习速率设为0.5, 如果调整为1,你会看到不同的结果
weight = -5 # 权重
bias = 3 # 偏置
# 计算一下初始权重和偏置值所带来的损失
print ('当前损失:',loss_function(X_train, y_train, weight, bias))
当前损失: 1.343795534906634
# 绘制当前的函数模型
plt.plot(X_train, y_train,'r.', label='Training data') # 显示训练集散点图
line_X = np.linspace(X_train.min(), X_train.max(), 500) # X值域
line_y = [weight*xx + bias for xx in line_X] # 假设函数y_hat
plt.plot(line_X,line_y,'b--', label='Current hypothesis' ) #显示当前拟合
plt.xlabel('Wechat Ads') # X轴Label
plt.ylabel('Sales') # y轴Label
plt.legend() # 显示图例
plt.show() # 显示绘图
# 根据初始参数值,进行梯度下降,也就是开始训练机器,拟合函数
loss_history, weight_history, bias_history = gradient_descent(
X_train, y_train, weight, bias, alpha, iterations)
plt.plot(loss_history,'g--',label='Loss Curve')
plt.xlabel('Iterations') # x轴Label
plt.ylabel('Loss') # y轴Label
plt.legend() # 显示图例
plt.show() # 显示损失曲线
# 绘制当前的函数模型
plt.plot(X_train, y_train,'r.', label='Training data') # 显示训练集散点图
line_X = np.linspace(X_train.min(), X_train.max(), 500) # X值域
# 关于weight_history[-1],这里的索引[-1],就代表迭代500次后的最后一个W值
line_y = [weight_history[-1]*xx + bias_history[-1] for xx in line_X] # 假设函数
plt.plot(line_X,line_y,'b--', label='Current hypothesis' ) # 显示当前函数
plt.xlabel('Wechat Ads') # x轴Label
plt.ylabel('Sales') # y轴Label
plt.legend() # 显示图例
plt.show() # 显示函数图像
print ('当前损失:',loss_function(X_train, y_train,
weight_history[-1], bias_history[-1]))
print ('当前权重:',weight_history[-1])
print ('当前偏置:',bias_history[-1])
当前损失: 0.006157577974318896
当前权重: 0.4721302015868674
当前偏置: 0.2707186935933173
print ('测试集损失:',loss_function(X_test, y_test,
weight_history[-1], bias_history[-1]))
测试集损失: 0.007776406285658269
# 同时绘制训练集和测试集损失曲线
loss_test ,a , b = gradient_descent(X_test, y_test, weight, bias, alpha, iterations)
plt.plot(loss_history,'g--',label='Traning Loss Curve')
plt.plot(loss_test,'r',label='Test Loss Curve')
plt.xlabel('Iterations') # x轴Label
plt.ylabel('Loss') # y轴Label
plt.legend() # 显示图例
plt.show()
# 设计Contour Plot动画
import matplotlib.animation as animation
# 初始化参数网格
theta0_vals = np.linspace(-2, 3, 100) # bias 范围
theta1_vals = np.linspace(-3, 3, 100) # weight 范围
J_vals = np.zeros((theta0_vals.size, theta1_vals.size))
# 计算损失函数网格
for t1, bias in enumerate(theta0_vals):
for t2, weight in enumerate(theta1_vals):
J_vals[t1, t2] = loss_function(X_train, y_train, weight, bias)
J_vals = J_vals.T
A, B = np.meshgrid(theta0_vals, theta1_vals)
C = J_vals
# 创建画布
fig = plt.figure(figsize=(12,5))
# 左图:数据拟合过程
ax1 = plt.subplot(121)
ax1.plot(X_train, y_train, 'ro', label='Training data')
ax1.set_title('Sales Prediction')
ax1.set_xlim(X_train.min()-X_train.std(), X_train.max()+X_train.std())
ax1.set_ylim(y_train.min()-y_train.std(), y_train.max()+y_train.std())
ax1.grid(axis='both')
ax1.set_xlabel("WeChat Ads Volumn (X1)")
ax1.set_ylabel("Sales Volumn (Y)")
ax1.legend(loc='lower right')
line, = ax1.plot([], [], 'b-', label='Current Hypothesis')
annotation = ax1.text(-2, 3, '', fontsize=12, color='green')
annotation.set_animated(True)
# 右图:损失函数等高线
ax2 = plt.subplot(122)
cp = ax2.contourf(A, B, C, levels=20, cmap='viridis')
plt.colorbar(cp, ax=ax2)
ax2.set_title('Loss Function Contour')
ax2.set_xlabel('Bias (θ₀)')
ax2.set_ylabel('Weight (θ₁)')
track, = ax2.plot([], [], 'r-', lw=1, alpha=0.5)
point, = ax2.plot([], [], 'ro', markersize=5)
plt.tight_layout()
# 动画初始化函数
def init():
line.set_data([], [])
track.set_data([], [])
point.set_data([], [])
annotation.set_text('')
return line, track, point, annotation
# 动画更新函数
def animate(i):
# 更新拟合直线
fit1_X = np.linspace(X_train.min()-X_train.std(), X_train.max()+X_train.std(), 1000)
fit1_y = bias_history[i] + weight_history[i] * fit1_X
line.set_data(fit1_X, fit1_y)
# 更新参数轨迹
track.set_data(bias_history[:i], weight_history[:i])
# 更新当前参数点
point.set_data([bias_history[i]], [weight_history[i]])
# 更新损失值显示
annotation.set_text(f'Iter: {i}\nCost = {loss_history[i]:.4f}')
return line, track, point, annotation
# 生成动画
anim = animation.FuncAnimation(
fig, animate, init_func=init,
frames=len(weight_history), # 确保与训练步数一致
interval=50, blit=True
)
# 保存动画
anim.save('linear_regression_training.gif', writer='pillow', fps=15)
# 显示Contour Plot动画
import io
import base64
from IPython.display import HTML
filename = 'linear_regression_training.gif'
video = io.open(filename, 'r+b').read()
encoded = base64.b64encode(video)
HTML(data='''<img src="data:image/gif;base64,{0}" type="gif" />'''.format(encoded.decode('ascii')))