回归算法是典型的监督学习模型之一。回归是一种统计学方法,用于根据样本数据
(
x
i
,
y
i
)
(\boldsymbol{x}_i,y_i)
(xi,yi),
i
=
1
,
2
,
⋯
,
m
i=1,2,\cdots,m
i=1,2,⋯,m,探究变量
x
\boldsymbol{x}
x与
y
y
y之间的关系。具体而言,回归模型的任务是找出拟合函数
F
(
x
)
F(\boldsymbol{x})
F(x),使得
y
i
≈
F
(
x
i
)
,
i
=
1
,
2
,
⋯
,
m
y_i\approx F(\boldsymbol{x}_i),i=1,2,\cdots,m
yi≈F(xi),i=1,2,⋯,m
并用
F
(
x
)
F(\boldsymbol{x})
F(x)来对新的输入
x
∈
R
n
\boldsymbol{x}\in\text{R}^n
x∈Rn预测对应的输出
y
∈
R
y\in\text{R}
y∈R。寻求拟合函数的过程,通常是选择一个具有待定参数
w
\boldsymbol{w}
w的函数
F
(
w
;
x
)
F(\boldsymbol{w};\boldsymbol{x})
F(w;x),其中
w
∈
R
p
,
p
∈
N
\boldsymbol{w}\in\text{R}^p,p\in\text{N}
w∈Rp,p∈N。然后计算使得
y
i
≈
F
(
w
0
;
x
i
)
,
i
=
1
,
2
,
⋯
,
m
y_i\approx F(\boldsymbol{w}_0;\boldsymbol{x}_i),i=1,2,\cdots,m
yi≈F(w0;xi),i=1,2,⋯,m最“合适”的参数
w
0
\boldsymbol{w}_0
w0——称为该学习模型的模式,并将
F
(
w
0
;
x
)
F(\boldsymbol{w}_0;\boldsymbol{x})
F(w0;x)作为拟合函数。
此处所谓最“合适”的模式,可以从不同的角度去认知。如果从欧氏空间的集合观点,可得出回归过程最常用的最小二乘法。给定序列
(
x
i
,
y
i
)
(\boldsymbol{x}_i,y_i)
(xi,yi),
i
=
1
,
2
,
⋯
,
m
i=1,2,\cdots,m
i=1,2,⋯,m,最小二乘法对所选含有待定参数
w
∈
R
p
\boldsymbol{w}\in\text{R}^p
w∈Rp的函数
F
(
w
;
x
)
F(\boldsymbol{w};\boldsymbol{x})
F(w;x),记
F
(
w
)
=
(
F
(
w
;
x
1
)
F
(
w
;
x
2
)
⋮
F
(
w
;
x
m
)
)
\boldsymbol{F}(\boldsymbol{w})=\begin{pmatrix} F(\boldsymbol{w};\boldsymbol{x}_1)\\F(\boldsymbol{w};\boldsymbol{x}_2)\\\vdots\\F(\boldsymbol{w};\boldsymbol{x}_m)\end{pmatrix}
F(w)=
F(w;x1)F(w;x2)⋮F(w;xm)
,
y
=
(
y
1
y
2
⋮
y
m
)
\boldsymbol{y}=\begin{pmatrix} y_1\\y_2\\\vdots\\y_m \end{pmatrix}
y=
y1y2⋮ym
,解无约束最优化问题
{
min
∥
F
(
w
)
−
y
∥
2
s.t
w
∈
R
p
,
\begin{cases} \min\quad\lVert\boldsymbol{F}(\boldsymbol{w})-\boldsymbol{y}\rVert^2\\ \text{s.t}\quad\boldsymbol{w}\in\text{R}^p \end{cases},
{min∥F(w)−y∥2s.tw∈Rp,
设
w
0
=
arg
min
w
∈
R
p
∥
F
(
w
)
−
y
∥
2
\boldsymbol{w}_0=\arg\min\limits_{\boldsymbol{w}\in\text{R}^p}\lVert\boldsymbol{F}(\boldsymbol{w})-\boldsymbol{y}\rVert^2
w0=argw∈Rpmin∥F(w)−y∥2,则
F
(
w
0
,
x
)
F(\boldsymbol{w}_0,\boldsymbol{x})
F(w0,x)即为所求的拟合函数。
需要指出的是,在回归模型中,无论是训练部分还是预测部分,都需要对数据作一些规范化处理。首先,将样本特征向量
x
i
=
(
x
i
1
x
i
2
⋮
x
i
n
)
∈
R
n
,
i
=
1
,
2
,
⋯
,
m
\boldsymbol{x}_i=\begin{pmatrix} x_{i1}\\x_{i2}\\\vdots\\x_{in} \end{pmatrix}\in\text{R}^n,i=1,2,\cdots,m
xi=
xi1xi2⋮xin
∈Rn,i=1,2,⋯,m将其组织成一个
m
×
n
m\times n
m×n矩阵
X
=
(
x
1
⊤
x
2
⊤
⋮
x
m
⊤
)
=
(
x
11
x
12
⋯
x
1
n
x
21
x
22
⋯
x
2
n
⋮
⋮
⋱
⋮
x
m
1
x
m
2
⋯
x
m
n
)
\boldsymbol{X}=\begin{pmatrix} \boldsymbol{x}_1^\top\\\boldsymbol{x}_2^\top\\\vdots\\\boldsymbol{x}_m^\top \end{pmatrix}=\begin{pmatrix} x_{11}&x_{12}&\cdots&x_{1n}\\ x_{21}&x_{22}&\cdots&x_{2n}\\ \vdots&\vdots&\ddots&\vdots\\ x_{m1}&x_{m2}&\cdots&x_{mn} \end{pmatrix}
X=
x1⊤x2⊤⋮xm⊤
=
x11x21⋮xm1x12x22⋮xm2⋯⋯⋱⋯x1nx2n⋮xmn
使得表达式展开时更简洁。其次,要对数据进行“标准化”处理。以消除不同量纲单位带来的数据偏差,并使数据指标处于同一数量级,更适合综合对比评价。此处以归一化作为规范化处理方式:
(1)训练阶段。此时,样本特征数据
x
i
,
i
=
1
,
2
,
⋯
,
m
\boldsymbol{x}_i,i=1,2,\cdots,m
xi,i=1,2,⋯,m且
m
>
1
m>1
m>1。对每个
1
≤
j
≤
n
1\leq j\leq n
1≤j≤n,第
j
j
j列数据
(
x
1
j
x
2
j
⋮
x
m
j
)
\begin{pmatrix} x_{1j}\\x_{2j}\\\vdots\\x_{mj} \end{pmatrix}
x1jx2j⋮xmj
表示
m
m
m个样本第
j
j
j个特征数据。计算每一列的最小值与最大值
min
x
j
=
min
i
{
x
i
j
}
\min x_j=\min\limits_{i}\{x_{ij}\}
minxj=imin{xij}及
max
x
j
=
max
i
{
x
i
j
}
\max x_j=\max\limits_{i}\{x_{ij}\}
maxxj=imax{xij}。并记
{
min
x
=
(
min
x
1
,
min
x
2
,
⋯
,
min
x
n
)
max
x
=
(
max
x
1
,
max
x
2
,
⋯
,
max
x
n
)
\begin{cases} \min\boldsymbol{x}=(\min x_1,\min x_2,\cdots,\min x_n)\\ \max\boldsymbol{x}=(\max x_1,\max x_2,\cdots,\max x_n) \end{cases}
{minx=(minx1,minx2,⋯,minxn)maxx=(maxx1,maxx2,⋯,maxxn)
令
Δ
x
=
max
x
−
min
x
=
(
max
x
1
−
min
x
1
,
max
x
2
−
min
x
2
,
⋯
,
max
x
n
−
min
x
n
)
\Delta\boldsymbol{x}=\max\boldsymbol{x}-\min\boldsymbol{x}=(\max x_1-\min x_1,\max x_2-\min x_2,\cdots,\max x_n-\min x_n)
Δx=maxx−minx=(maxx1−minx1,maxx2−minx2,⋯,maxxn−minxn),以
x
i
⊤
−
min
x
Δ
x
=
(
x
1
−
min
x
1
max
x
1
−
min
x
1
,
x
2
−
min
x
2
max
x
2
−
min
x
2
,
⋯
,
x
n
−
min
x
n
max
x
n
−
min
x
n
)
\frac{\boldsymbol{x}_i^\top-\min\boldsymbol{x}}{\Delta\boldsymbol{x}}=\left(\frac{x_1-\min x_1}{\max x_1-\min x_1},\frac{x_2-\min x_2}{\max x_2-\min x_2},\cdots,\frac{x_n-\min x_n}{\max x_n-\min x_n}\right)
Δxxi⊤−minx=(maxx1−minx1x1−minx1,maxx2−minx2x2−minx2,⋯,maxxn−minxnxn−minxn)
作为
x
i
⊤
\boldsymbol{x}_i^\top
xi⊤归一化后的向量,
i
=
1
,
2
,
⋯
,
m
i=1,2,\cdots,m
i=1,2,⋯,m。
相仿地,对接收到的标签数据
y
i
y_i
yi,
i
=
1
,
2
,
⋯
,
m
i=1,2,\cdots,m
i=1,2,⋯,m,记
min
y
=
min
i
{
y
i
}
,
max
y
=
max
i
{
y
i
}
\min y=\min\limits_{i}\{y_i\},\max y=\max\limits_{i}\{y_i\}
miny=imin{yi},maxy=imax{yi},以
y
i
−
min
y
max
y
−
min
y
\frac{y_i-\min y}{\max y-\min y}
maxy−minyyi−miny
作为
y
i
y_i
yi归一化后的值。数据归一化后,其值均介于0,1之间。
(2)预测阶段。此时,只有新样本的特征向量
x
i
\boldsymbol{x}_i
xi,
i
=
1
,
2
,
⋯
,
m
i=1,2,\cdots,m
i=1,2,⋯,m。进行归一化时,需要用训练时算得的
min
j
\min_j
minj和
max
j
\max_j
maxj进行计算。这是因为,我们认为训练数据是总体的简单样本,其统计特征值(最小值、最大值)表示总体的近似分布。而检测数据也来自同一总体,出于一致性考虑,使用训练时取得的最大值、最小值对检测数据作归一化计算。
回归模型中,若拟合函数选择为线性函数,即
y
=
F
(
w
;
x
)
=
∑
i
=
1
n
w
i
x
i
+
w
n
+
1
=
(
x
⊤
,
1
)
w
y=F(\boldsymbol{w};\boldsymbol{x})=\sum_{i=1}^nw_ix_i+w_{n+1}=(\boldsymbol{x}^\top,1)\boldsymbol{w}
y=F(w;x)=i=1∑nwixi+wn+1=(x⊤,1)w
则称为线性回归模型。其中,
w
1
,
w
2
,
⋯
,
w
n
w_1,w_2,\cdots,w_n
w1,w2,⋯,wn为
x
1
,
x
2
,
⋯
,
x
n
x_1,x_2,\cdots,x_n
x1,x2,⋯,xn的加权和
∑
i
=
1
n
w
i
x
i
\sum\limits_{i=1}^nw_ix_i
i=1∑nwixi的系数,
w
n
+
1
w_{n+1}
wn+1为偏移量。线性回归模型可图示化为
下面,将线性回归模型定义成类LineModel。
import numpy as np #导入numpy
from scipy.optimize import minimize #导入minimize
class LineModel(): #线性模型
def xnormalize(self, x, trained): #样本特征数据归一化方法
if not trained: #训练前
xmin = np.min(x,axis = 0) #按列计算最小值
xmax = np.max(x,axis = 0) #按列计算最大值
self.xmin = xmin #记录计算结果
self.xmax = xmax
else: #测试或预测前
xmin = self.xmin #使用训练时记录的数据
xmax = self.xmax
return (x - xmin) / (xmax - xmin)
def ynormalize(self, y, trained): #样本标签归一化方法
if not trained: #训练前
self.ymin = np.min(y) #记录最小值
self.ymax = np.max(y) #记录最大值
return (y - self.ymin)/(self.ymax - self.ymin)
def pretreat(self, x, y = None, trained = False): #数据预处理函数
if isinstance(y, np.ndarray): #需处理样本标签
y = self.ynormalize(y, trained) #归一化标签
if not isinstance(x, np.ndarray): #一个一元样本
x = np.array([x]).reshape(1, 1)
else:
if len(x.shape) == 1: #一维数组
if self.scalar: #多个一元样本
x = x.reshape(x.size,1)
else: #一个多元样本
x = x.reshape(1,x.size)
x = self.xnormalize(x, trained) #归一化样本特征
m = x.shape[0] #样本个数
return np.hstack((x, np.ones(m).reshape(m, 1))), y #x添加一列1
def fit(self, X, Y, w=None): #训练方法
print("训练中...,稍候")
self.scalar = (len(X.shape)==1) #记录一元样本特征
self.A, self.y = self.pretreat(X, Y) #预处理训练数据
p = self.patternlen() #计算模式长度
if w == None: #未传递初始参数向量w
w = np.random.random(p) #随机产生初始向量w
else:
if type(w) == int or type(w) == float: #常数
w=np.array([w]*p)
res = minimize(self.obj, w, #解最优化问题
method = conjFR,options={'gtol':1e-4})
self.pattern = res.x #记录模式
print("%d次迭代后完成训练。"%res.nit)
def predict(self, X): #预测方法
X, _ = self.pretreat(X, trained=True) #归一化样本特征数据
yp = self.F(self.pattern, X) #计算拟合函数值
if yp.size == 1: #单样本预测值
yp = yp[0]
return yp*(self.ymax - self.ymin) + self.ymin
def obj(self,w): #目标函数
return np.linalg.norm(self.F(w, self.A) - (self.y))**2
def patternlen(self): #模式长度函数
return self.A.shape[1]
def F(self, w, x): #线性拟合函数
return np.matmul(x,w)
借助代码内的注释信息,看官不难理解程序。需要注意的是,训练函数fit中第43~44行调用scipy.optimize的minimize函数(第2行导入)计算目标函数obj的最优解,所用的方法conjFR是我们在博文《最优化方法Python计算:非二次型共轭梯度算法》中定义的。
例1 为研究某一化学反应过程中,温度
x
x
x(℃)对产品得率
y
y
y(%)的影响,测得数据如下
x x x | 100 100 100 | 110 110 110 | 120 120 120 | 130 130 130 | 140 140 140 | 150 150 150 | 160 160 160 | 170 170 170 | 180 180 180 | 200 200 200 |
---|---|---|---|---|---|---|---|---|---|---|
y y y | 45 45 45 | 51 51 51 | 54 54 54 | 61 61 61 | 66 66 66 | 70 70 70 | 74 74 74 | 78 78 78 | $85 | 89 89 89 |
记
x
=
(
100
,
110
,
120
,
130
,
140
,
150
,
160
,
170
,
180
,
200
)
⊤
\boldsymbol{x}=(100,110,120,130,140,150,160,170,180,200)^\top
x=(100,110,120,130,140,150,160,170,180,200)⊤,
y
=
(
45
,
51
,
54
,
61
,
66
,
70
,
74
,
78
,
85
,
89
)
⊤
\boldsymbol{y}=(45,51,54,61,66,70,74,78,85,89)^\top
y=(45,51,54,61,66,70,74,78,85,89)⊤。试用上述程序中定义的LinearModel类创建一个线性回归模型,并用数据
x
\boldsymbol{x}
x和
y
\boldsymbol{y}
y训练模型,预测温度
x
=
147
x=147
x=147时的得率。
解:下列代码完成计算
import numpy as np #导入numpy
x=np.array([100, 110, 120, 130, 140, 150, 160, 170, 180, 200]) #样本特征数据
y=np.array([45, 51, 54, 61, 66, 70, 74, 78, 85, 89]) #标签数据
chemical = LineModel() #构造线性回归模型
chemical.fit(x, y) #训练
pattern = chemical.pattern #最优模式
print('模式:%s'%pattern)
x1 = 147 #新的样本特征
print('对温度x1=%d,预测得率y=%.2f'%(x1,chemical.predict(x1)) + '%')
运行程序,输出
训练中...,稍候
3次迭代后完成训练。
模式:[1.03354978 0.03138529]
对温度x1=147,预测得率y=67.75%
写博不易,敬请支持:
如果阅读本文于您有所获,敬请点赞、评论、收藏,谢谢大家的支持!