一、线性回归前期准备
1. 什么是线性回归?
当我们求解一个问题,这个问题跟n个因素相关。我们需要求出每一个因素的权重
θ
i
(
1
<
=
i
<
=
n
)
\theta_{i}(1<=i<=n)
θi(1<=i<=n) ,这个过程就叫做线性回归。
比如我们求一个地区的房价,最起码和房价相关的因素应该有 地区、面积、户型 等,我们就可以得到下面表达式。其中
θ
0
\theta_{0}
θ0 是常数项的一个权值,也就是我们初中的y=ax+b中的b,而
x
0
x_{0}
x0 只是为了配合前面的b 而特地添加的,并无实际的意义,一般设为1,这样 1*b=b。不影响最后的结果。而求解
θ
i
\theta_{i}
θi的过程我们就叫做线性回归。
Y p r i c e = θ 0 x 0 + θ 1 x 1 + θ 2 x 2 + . . . + θ n x n Y_{price}=\theta_{0}x_{0}+\theta_{1}x_{1}+\theta_{2}x_{2}+...+\theta_{n}x_{n} Yprice=θ0x0+θ1x1+θ2x2+...+θnxn
2. 什么是代价函数?
在机器学习中,没有完美的权值。不管我们随机生成的初始
θ
i
\theta_i
θi值和最终通过算法迭代求的最终值都有误差。而代价函数就是为了去表示这个误差。一般的代价函数表示如下。其中m是值数据的条数。y是指真实的结果。所以cost函数可以表示为:
c
o
s
t
=
1
/
m
∑
i
=
1
m
(
y
i
−
h
θ
(
x
i
)
)
2
cost=1/m\sum_{i=1}^m(y_{i}-h_{\theta}(x_{i}))^{2}
cost=1/mi=1∑m(yi−hθ(xi))2
其中:
h
θ
(
x
)
=
θ
0
x
0
+
θ
1
x
1
+
θ
2
x
2
+
.
.
.
+
θ
n
x
n
h_{\theta}(x)=\theta_{0}x_{0}+\theta_{1}x_{1}+\theta_{2}x_{2}+...+\theta_{n}x_{n}
hθ(x)=θ0x0+θ1x1+θ2x2+...+θnxn
Q:为什么
(
y
i
−
h
θ
(
x
i
)
)
2
(y_{i}-h_{\theta}(x_{i}))^{2}
(yi−hθ(xi))2为什么要平方?
A:因为(1)我们通过权值求出来的值和真实值大小不好确定,所以我们采取平方的方法来减少因为正负带来的误差抵消。(2)更为重要的是,这个函数是一个凸函数,在函数的取值范围内处处可导,便于后面求梯度。
Q:为什么要1/m?
A:我们最终求得是一个平均代价。其中这里面并没有太多的关系,因为我们的目标是使cost变小。我们给cost乘一个常数并不会影响最后的变大变小的趋势,不除m或者乘100也是没有关系的。
3. 什么是偏导?
偏导就是我们一个方程有多个未知数,我们不可能对所有的未知数全部求导。我们一次只能对一个参数进行求导。例:
Y
(
θ
1
,
θ
2
)
=
m
∗
θ
1
+
b
∗
θ
2
d
Y
d
θ
1
=
m
Y(\theta_{1},\theta_{2})=m*\theta_{1}+b*\theta_{2}\\ \frac{dY}{d\theta_{1}}=m
Y(θ1,θ2)=m∗θ1+b∗θ2dθ1dY=m
4. 什么是梯度下降?
梯度下降是运用了偏微分的思想来求解 θ i \theta_{i} θi的一种思路。在多变量函数中,梯度是一个向量,向量有方向,梯度的方向就指出了函数在给定点的上升最快的方向,而梯度取负号就是整个函数下降最快的方向。这样我们就可以利用梯度下降的方式通过改变参数的方式让我们的代价函数(真实值和我们求出来的值之差)逐渐缩小。这样我们就求出了最优的 θ i \theta_i θi。
二、梯度下降
1. 为什么要用梯度下降?
我们前文讲到代价函数,我们的目标是使这个代价函数最小(趋近于0).在最小的时候,我们就可以拿到此时的
θ
\theta
θ, 这就是我们想要的权值啊。但是怎么让一个多元的函数最小呢。
1.暴力(遍历所有的可能)不可取,因为数据维度(有很多个相关参数)太大,参数取值无范围。这样岂不是要暴力到猴年马月。
2.智取。我们可以想象一元两次方程
y
=
x
2
+
4
y=x^2+4
y=x2+4我们可以先求导数,然后令导数等于0,这样我们就可以在那一点求得驻点,然后我们可以分析这个点是不是最大值或者最小值。
我们也可以在某一点求导数,比如我们在x=2的地方求导=4 我们可以理解为在2后面的斜率是4.我们想让我们的y小,我们就去往反方向走。比如我们在-2求得要往正方向走才能使y变小。
推广到多元变量,多元函数中要想让代价函数值变得更小。我们也应该沿着每一个参数的负梯度的方向走。
2.怎么运用梯度下降的方法来求最终的权值
- 首先我们先随机选取一组初始权值来求最后的结果
- 循环迭代,每次我们往梯度下降的方向改变一点点。
梳理公式
h
θ
(
x
i
)
=
θ
0
x
0
+
θ
1
x
1
+
θ
2
x
2
+
.
.
.
+
θ
n
x
n
(拟合函数)
h_{\theta}(x_{i})=\theta_{0}x_{0}+\theta_{1}x_{1}+\theta_{2}x_{2}+...+\theta_{n}x_{n} \tag{拟合函数}
hθ(xi)=θ0x0+θ1x1+θ2x2+...+θnxn(拟合函数)
J
(
θ
)
=
1
/
m
∑
i
=
1
m
(
y
i
−
h
θ
(
x
i
)
)
2
(代价函数)
J(\theta)=1/m\sum_{i=1}^m(y^{i}-h_{\theta}(x^{i}))^{2} \tag{代价函数}
J(θ)=1/mi=1∑m(yi−hθ(xi))2(代价函数)
公式说明:在拟合公式中,
θ
i
\theta_{i}
θi当前已知的,因为已经假设了或者处在迭代期间,不管正确与否,反正当前已知。而x是要带入公式去求最后的拟合值的。(其中x也是已知的,因为这是我们的数据啊,因为这里我们要带入每一份数据的x)。我们姑且在这里面称x是变量。
在代价函数中,是关于
θ
\theta
θ的一个多元函数。因为在这里面我们最终就是要求它。
推导过程
目标:我们要使代价函数尽可能的小
方法:我们采取对每一个权值往负梯度方向走。
- 我们首先为代价函数的每一个
θ
\theta
θ求偏导
δ J ( θ ) δ θ k = 1 m ∑ i = 1 m ( y i − h θ ( x i ) ) x k i (通式) \frac {\delta J(\theta)}{\delta \theta_k}=\frac{1}{m}\sum_{i=1}^m(y^{i}-h_{\theta}(x^{i}))x_k^i \tag{通式} δθkδJ(θ)=m1i=1∑m(yi−hθ(xi))xki(通式) δ J ( θ ) δ θ 1 = 1 m ∑ i = 1 m ( y i − h θ ( x i ) ) x 1 i ( θ 1 ) \frac {\delta J(\theta)}{\delta \theta_1}=\frac{1}{m}\sum_{i=1}^m(y^{i}-h_{\theta}(x^{i}))x_1^i \tag{$\theta_1$} δθ1δJ(θ)=m1i=1∑m(yi−hθ(xi))x1i(θ1)
δ J ( θ ) δ θ = < δ J ( θ ) δ θ 0 , δ J ( θ ) δ θ 1 , δ J ( θ ) δ θ 2 . . . > (整体求出来是n+1个向量) \frac {\delta J(\theta)}{\delta \theta}=<\frac {\delta J(\theta)}{\delta \theta_0},\frac {\delta J(\theta)}{\delta \theta_1},\frac {\delta J(\theta)}{\delta \theta_2}...>\tag{整体求出来是n+1个向量} δθδJ(θ)=<δθ0δJ(θ),δθ1δJ(θ),δθ2δJ(θ)...>(整体求出来是n+1个向量)
2.正式推导(t+1是t的下一轮)
J ( θ t + 1 ) − J ( θ t ) = ( θ t + 1 − θ t ) J ( θ t ) ′ + R ( i g n o r e ) (泰勒公式) J(\theta_{t+1})-J(\theta_{t})=(\theta_{t+1}-\theta_t)J(\theta_t)'+R(ignore)\tag{泰勒公式} J(θt+1)−J(θt)=(θt+1−θt)J(θt)′+R(ignore)(泰勒公式)
若想
J ( θ t + 1 ) < J ( θ t ) J(\theta_{t+1})<J(\theta_{t}) J(θt+1)<J(θt)
则
( θ t + 1 − θ t ) J ( θ t ) ′ < 0 (\theta_{t+1}-\theta_t)J(\theta_t)'<0 (θt+1−θt)J(θt)′<0
可令
( θ t + 1 − θ t ) = − J ( θ t ) ′ (\theta_{t+1}-\theta_t)=-J(\theta_t)' (θt+1−θt)=−J(θt)′
可得
J ( θ t + 1 ) − J ( θ t ) = − ( J ( θ t ) ′ ) 2 < = 0 J(\theta_{t+1})-J(\theta_{t})=-(J(\theta_t)')^2<=0 J(θt+1)−J(θt)=−(J(θt)′)2<=0
所以在 J ( θ t + 1 ) < J ( θ t ) J(\theta_{t+1})<J(\theta_{t}) J(θt+1)<J(θt)的目标下我们得到了 ( θ t + 1 − θ t ) = − J ( θ t ) ′ (\theta_{t+1}-\theta_t)=-J(\theta_t)' (θt+1−θt)=−J(θt)′。即 θ t + 1 = θ t − J ( θ t ) ′ \theta_{t+1}=\theta_t-J(\theta_t)' θt+1=θt−J(θt)′
由于泰勒公式只在 θ t + 1 − θ t \theta_{t+1}-\theta_{t} θt+1−θt很小时才成立。我们不用刻意去让他们的差趋近于0,那样运算时间太长了。我们只要小一点就可以。所以我们引入 α \alpha α学习率(步长)来让 θ \theta θ变得慢。在机器学习中, α \alpha α通常取0.001,0.0001.
我们得到了我们想要的使权值变好的公式。
Q:既然我们目标是使代价函数变成0,为什么不直接求出来最优解,而是一步步缓慢的去迭代。
A:(我的理解)这就变成了一个多元的函数求解的问题,我们可能需要带入一定多的数据进行求解。但这些数据中难免有噪音数据。即使我们用全部的数据都去求 然后最后实现类似平均值之类的话,也需要同时大量的内存消耗,还有就是求出来的解并不一定是全局的最优解,有可能是一个局部的最优解。欢迎讨论。
代码实现。
代价函数
代价函数是为了计算总的一个损失的。是为了画图用。
def computerCost(X, y, theta):
# X是一个m行、n列的矩阵。y是一个m行1列的矩阵,theta是 1 行,n 列的矩阵。
# 计算代价,主要用于画图时模拟一下。
# 实现对应元素相乘,有2种方式,一个是np.multiply(),另外一个是 *
return (np.transpose(X * theta - y)) * (X * theta - y) / (2 * len(y))
线性回归
def linearRegression(data, alpha, num_iters, isUnification):
X = data[:, 0:-1] # X是前n-1列,即自变量的列
y = data[:, -1] # y对应最后一列,即函数值
m = len(y) # 总的数据条数
col = data.shape[1] # data的列数
if isUnification:
X, mu, sigma = featureNormaliza(X) # 归一化 均值mu,标准差sigma,
if X.shape[1]<=2 and X.shape[1]>=1:
plot_data(data) # 绘制数据的归一化效果
X = np.hstack((np.ones((m, 1)), X)) # 在X前加一列1,这个是为了表示后面那个常数的,*1还是等于原来的数
theta = np.zeros((col, 1)) # 初始化权值向量
y = y.reshape(-1, 1) # 将行向量转化为列
theta, history = gradientDescent(X, y, theta, alpha, num_iters)
plotJ(history, num_iters)
return theta # 返回和学习的结果theta
梯度下降
def gradientDescent(X, y, theta, alpha, num_iters):
# X是一个m行3列的一个矩阵、y是每一个Xi的函数值、theta是权值、alpha是每一次的步长(学习率)、num是迭代的次数
m = len(y)
n = len(theta)
history = np.zeros((num_iters, 1)) # 记录每次迭代计算的代价值
for i in range(num_iters): # 遍历迭代次数
h = np.dot(X, theta) # 计算内积,matrix可以直接乘
a = ((alpha / m) * (np.dot(np.transpose(X), h - y)))
theta = np.mat(theta - a) # 梯度的计算
history[i] = computerCost(X, y, theta) # 调用计算代价函数
if (i % 100000 == 0):
print('总共', num_iters, ' 当前迭代了', i / num_iters * 100, '%了')
return theta, history
总代码
import numpy as np
from matplotlib import pyplot as plt
from matplotlib.font_manager import FontProperties
font = FontProperties(fname=r"c:\\windows\\fonts\\simsun.ttc", size=14) # 解决windows环境下画图汉字乱码问题
X = np.random.normal(0, 5, size=(100, 4))
theta = np.array([3, 4, 5, 6])
b = 5
y = np.array([np.sum(theta * i) + b + np.random.normal(4, 10) for i in X])
# 画每次迭代代价的变化图
def plotJ(history, num_iters):
x = np.arange(1, num_iters + 1)
plt.plot(x, history)
plt.xlabel(u"迭代次数", fontproperties=font) # 注意指定字体,要不然出现乱码问题
plt.ylabel(u"代价值", fontproperties=font)
plt.title(u"代价随迭代次数的变化", fontproperties=font)
plt.show()
# 归一化feature
def featureNormaliza(X):
X_norm = np.array(X) # 将X转化为numpy数组对象,才可以进行矩阵的运算
# 定义所需变量
mu = np.mean(X_norm, 0) # 求每一列的平均值(0指定为列,1代表行)
sigma = np.std(X_norm, 0) # 求每一列的标准差
for i in range(X.shape[1]): # 遍历列
X_norm[:, i] = (X_norm[:, i] - mu[i]) / sigma[i] # 归一化
return X_norm, mu, sigma
def computerCost(X, y, theta):
# X是一个m行、n列的矩阵。y是一个m行1列的矩阵,theta是 1 行,n 列的矩阵。
# 计算代价,主要用于画图时模拟一下。
# 实现对应元素相乘,有2种方式,一个是np.multiply(),另外一个是 *
a = (np.transpose(X * theta - y))
b = (X * theta - y)
c = a * b
return c / (2 * len(y))
# 梯度下降算法
def gradientDescent(X, y, theta, alpha, num_iters):
# X是一个m行3列的一个矩阵、y是每一个Xi的函数值、theta是权值、alpha是每一次的步长(学习率)、num是迭代的次数
m = len(y)
n = len(theta)
history = np.zeros((num_iters, 1)) # 记录每次迭代计算的代价值
for i in range(num_iters): # 遍历迭代次数
h = np.dot(X, theta) # 计算内积,matrix可以直接乘
a = ((alpha / m) * (np.dot(np.transpose(X), h - y)))
theta = np.mat(theta - a) # 梯度的计算
history[i] = computerCost(X, y, theta) # 调用计算代价函数
if (i % 100000 == 0):
print('总共', num_iters, ' 当前迭代了', i / num_iters * 100, '%了')
return theta, history
def linearRegression(X, y, alpha, num_iters, isUnification):
m = len(y) # 总的数据条数
col = X.shape[1] + 1 # data的列数
if isUnification:
X, mu, sigma = featureNormaliza(X) # 归一化 均值mu,标准差sigma,
X = np.hstack((np.ones((m, 1)), X)) # 在X前加一列1,这个是为了表示后面那个常数的,*1还是等于原来的数
theta = np.zeros((col, 1)) # 初始化权值向量
y = y.reshape(-1, 1) # 将行向量转化为列
theta, history = gradientDescent(X, y, theta, alpha, num_iters)
plotJ(history, num_iters)
return theta # 返回和学习的结果theta
# 测试学习效果(预测)
def predict(mu, sigma, theta): # 每一列的平均值,平均差,权值
predict = np.array([1650, 3])
norm_predict = (predict - mu) / sigma
final_predict = np.hstack((np.ones((1)), norm_predict))
result = np.dot(final_predict, theta) # 预测结果
return result
if __name__ == "__main__":
isUnification = False # 是否归一化, 归一化后 测试数据也要归一化才能进行预测。
alpha = 0.0001 # 学习率
num_iter = 5000 # 迭代次数
theta = linearRegression(X, y, alpha, num_iter, isUnification)
for i in range(len(theta)): # 输出权值
print('权值', i, '=', theta[i][0])
数据生成代码
GitHub
如果感觉还可以麻烦给个star。
问答
- Q:梯度下降?批量梯度下降?随机梯度下降?小批量梯度下降?
A:梯度下降就是一个算法总称,或者你可以理解成一个批量梯度下降的一个特殊形式,只有一条数据。批量梯度下降是我们本文重点讲解的。有很多条数据来求最终的权值。
随机梯度下降是批量梯度下降的一种缩减形式(一般用的最多)。
h = np.dot(X, theta) # 计算内积,matrix可以直接乘
a = ((alpha / m) * (np.dot(np.transpose(X), h - y)))
theta = np.mat(theta - a) # 梯度的计算
上面代码是批量梯度下降的,我们每次计算下一轮的权值都要把m条数据统统计算一遍,才能得到下一轮的权值。这样太费时间了,随机梯度下降就是我们只取出其中一条数据来计算下一轮的权值,因为只取一条,所以风险很大,因为这一条数据也许是一个激进的数据,这一轮的权值有可能对比上一轮的权值会差别很大。但是有风险肯定有优点,他说不准就会脱离一个局部最优解,进入全局最优解。
而小批量梯度下降是 批量和随机的一个结合。我们随机选10(小批量)条数据,不用选全部的,也没必要冒风险只选一条数据。
注意事项
- 我们什么时候用归一化?
归一化就是把数据 统一归到一个小范围内
当我们的数据不同维度的权重差别很大,比如,我们y=ax1+bx2+c 我们可能知道a这个权值很大,也就是说x1对y的影响很大。是b的1000倍,这时候我们就需要把数据进行归一化,因为 我们初始权值是一样的 我们的步长也是固定不变的,如果不进行归一化,每次对于b来说可能改变的太多了,而对于a来说只是一点点。
归一化之后,迭代次数要更多才能达到类似的效果,当数据值太大如果不归一下,可能就得不到最优的解。对比
没归一化1000次迭代就差不多了
有了归一化要迭代百万级才能(因为都缩小到很小的范围内,导数都很小,改变慢):
- 学习率怎么选
我们的数据值大小很大时(>100),我们需要把学习率设置的小一点 10E-6就可,否则一点点的改变也许就跳过了最优解进入了发散的状态。