1.Python代码
#!/usr/bin/env python3
# encoding: utf-8
'''
@file: larsWine.py
@time: 2020/6/13 0013 8:59
@author: Jack
@contact: jack18588951684@163.com
'''
import urllib.request
import numpy as np
from sklearn import datasets, linear_model
from math import sqrt
import matplotlib.pyplot as plt
## 读取数据集
target_url = ("http://archive.ics.uci.edu/ml/machine-learning-databases/wine-quality/winequality-red.csv")
data = urllib.request.urlopen(target_url)
xList = []
labels = []
names = []
firstLine = True
for line in data:
if firstLine:
names = str(line, encoding='utf-8').strip().split(";")
firstLine = False
else:
row = str(line, encoding='utf-8').strip().split(";")
labels.append(float(row[-1]))
row.pop()
floatRow = [float(num) for num in row]
xList.append(floatRow)
nrows = len(xList)
ncols = len(xList[0])
## 计算均值方差
xMean = []
xSD = []
for i in range(ncols):
col = [xList[j][i] for j in range(nrows)]
mean = sum(col) / nrows
xMean.append(mean)
colDiff = [(xList[j][i] - mean) for j in range(nrows)]
sumSq = sum([colDiff[i] * colDiff[i] for i in range(nrows)])
stdDev = sqrt(sumSq / nrows)
xSD.append(stdDev)
## 归一化属性和标签
xNormalized = []
for i in range(nrows):
rowNormalized = [(xList[i][j] - xMean[j]) / xSD[j] for j in range(ncols)]
xNormalized.append(rowNormalized)
meanLabel = sum(labels) / nrows
sdLabel = sqrt(sum([(labels[i] - meanLabel) * (labels[i] - meanLabel) for i in range(nrows)]) / nrows)
labelNormalized = [(labels[i] - meanLabel) / sdLabel for i in range(nrows)]
## 权重矩阵初始化
beta = [0.0] * ncols
betaMat = []
betaMat.append(list(beta))
nSteps = 350 ## 迭代次数
stepSize = 0.004 ## 步长
for i in range(nSteps):
residuals = [0.0] * nrows
for j in range(nrows):
## 预测值
labelsHat = sum([xNormalized[j][k] * beta[k] for k in range(ncols)])
## 误差
residuals[j] = labelNormalized[j] - labelsHat
## 计算相关度
corr = [0.0] * ncols
for j in range(ncols):
corr[j] = sum([xNormalized[k][j] * residuals[k] for k in range(nrows)]) / nrows
iStar = 0
corrStar = corr[0]
for j in range(1, (ncols)):
if abs(corrStar) < abs(corr[j]):
iStar = j
corrStar = corr[j]
## 权重更新
beta[iStar] += stepSize * corrStar / abs(corrStar)
betaMat.append(list(beta))
for i in range(ncols):
coefCurve = [betaMat[k][i] for k in range(nSteps)]
xaxis = range(nSteps)
plt.plot(xaxis, coefCurve)
plt.xlabel("Steps Taken")
plt.ylabel("Coefficient Values")
plt.show()
2.代码说明
代码主要包括3个部分:
(1)读入数据以及列名称,将数据转换为属性列表(Python的list对象,每个元素对应一个样例,该样例使用包含属性的list对象表示)以及标签。
(2)对属性以及标签进行归一化。
(3)对问题进行求解获得结果 β* β0 。
第一部分用于读入整个文件,将文件头分离出来,使用符号“;”获得文件头分隔的属性名列表,将剩下的行切分为浮点数列表,将所有记录封装为列表(一行记录)的列表,将标签封装为列表。对这些数据结构使用 Python 的列表(list)类型,是因为算法要对行和列进行遍历,而Pandas的数据帧类型对这个目的来说太慢。
第二部分属性归一化的目的是将属性值转换为同等尺度,从而可以方便地绘制以及充分填充坐标。在惩罚线性回归中,归一化一般作为第一步,原因大致相同。算法的每一步都会固定增加β变量值。属性尺度不同,固定增量对不同属性的影响也不同。同样,如果改变一个属性的尺度(从英里变为英尺),答案会显著不同。因为这些原因,惩罚线性回归一般归约到0均值(减去平均值)、单位方差(用标准方差去除)。
第 三部分求解 β、β0*。因为算法运行在归一化的变量上,因此不需要截距β0 。截距β0 一般用于表示标签值与加权属性值之间的差异。因为所有属性已经被归约到0均值,所以标签与加权属性的期望没有偏差,β0 * 就不再需要了。注意到有2个beta相关的列表被初始化。一个称作beta,其中元素个数与属性个数相同。另一个是类似于矩阵 - 列表的列表 - 用于存储算法每一步生成的beta列表。这些都涉及了惩罚线性回归以及现代机器学习算法的重要概念。
3.结果分析
基本迭代过程只包含几行代码,从 for 循环开始迭代 n 步。迭代的起始点是 β 的一个值。第一遍,所有值被设为 0。后续遍历,使用上一次迭代的结果。每次迭代包括 2 步。首先,使用 β 计算误差(指真实标签与预测结果的差异)。在该例中,预测结果值等于属性值乘以系数然后加和。第 2 步找到每个属性与误差的关联,从而决定哪个属性对降低误差贡献最大。2 个变量的关联值等于归一化(减去均值后除以标准差)值的乘积。如果 2 个变量完全相关,比如一个变量是另一个变量的缩放版,那以正的缩放对应于正关联,负的缩放对应于负关联。如果 2 个变量相互独立,那么他们的关联系数为 0。列表 corr 包含每个属性的关联计算结果。一旦关联计算完成,决定哪个属性与误差有最大关联(绝对值最大)就变得简单。β列表中的对应元素会增加一点。如果关联为正,增量就为正,否则增量为负。之后利用 β的新值再进行迭代计算。 算法生成的结果如图中的曲线所示。
对图进行理解可以按照如下的方法:
沿着迭代方向想象有一个点,在该点上,一条垂直线会穿过所有系数曲线。垂直线与系数曲线相交的值是算法在该步得到的系数。如果生成曲线需要 350 步,对应就会有350 个系数集合。注意到对于前 25 步,只有一个系数值非 0。这对应于套索回归(L1正则化)的稀疏属性。随着迭代步数的增加,第2个变量开始出现,这样的过程一直持续到所有变量被包含到解中。系数远离 0的过程可以看作是变量重要性的反映。在某些情况下,如果要去掉一些变量的预测,就应该尽量去掉后面出现的而不是开始出现的变量。
生成变量重要性排序是惩罚线性回归模型的一个重要特征。这使得该模型成为模型初步研究中的一个便捷的工具,因为它们将帮助构建关于保留哪个变量、丢弃哪个变量的决定,这个过程也称作特征工程。