1.python-梯度下降法 做重回归分析
准备好数据文件csv：
直接上代码：
import numpy as np
from matplotlib import pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
data=np.genfromtxt('./漫画数据2.csv',delimiter=',')
x_data=data[:,:-1]
y_data=data[:,2]
#定义学习率、斜率、截据
#设方程为y=theta1x1+theta2x2+theta0
lr=0.00001
theta0=0
theta1=0
theta2=0
#定义最大迭代次数，因为梯度下降法是在不断迭代更新k与b
epochs=10000
#定义最小二乘法函数-损失函数（代价函数）
def compute_error(theta0,theta1,theta2,x_data,y_data):
totalerror=0
for i in range(0,len(x_data)):#定义一共有多少样本点
totalerror=totalerror+(y_data[i]-(theta1*x_data[i,0]+theta2*x_data[i,1]+theta0))**2
return totalerror/float(len(x_data))/2
#梯度下降算法求解参数
def gradient_descent_runner(x_data,y_data,theta0,theta1,theta2,lr,epochs):
m=len(x_data)
for i in range(epochs):
theta0_grad=0
theta1_grad=0
theta2_grad=0
for j in range(0,m):
theta0_grad-=(1/m)*(-(theta1*x_data[j,0]+theta2*x_data[j,1]+theta2)+y_data[j])
theta1_grad-=(1/m)*x_data[j,0]*(-(theta1*x_data[j,0]+theta2*x_data[j,1]+theta0)+y_data[j])
theta2_grad-=(1/m)*x_data[j,1]*(-(theta1*x_data[j,0]+theta2*x_data[j,1]+theta0)+y_data[j])
theta0=theta0-lr*theta0_grad
theta1=theta1-lr*theta1_grad
theta2=theta2-lr*theta2_grad
return theta0,theta1,theta2
#进行迭代求解
theta0,theta1,theta2=gradient_descent_runner(x_data,y_data,theta0,theta1,theta2,lr,epochs)
print('结果：迭代次数：{0} 学习率：{1}之后 a0={2},a1={3},a2={4},代价函数为{5}'.format(epochs,lr,theta0,theta1,theta2,compute_error(theta0,theta1,theta2,x_data,y_data)))
print("多元线性回归方程为:y=",theta1,"X1+",theta2,"X2+",theta0)
#画图
ax=plt.figure().add_subplot(111,projection='3d')
ax.scatter(x_data[:,0],x_data[:,1],y_data,c='r',marker='o')
x0=x_data[:,0]
x1=x_data[:,1]
#生成网格矩阵
x0,x1=np.meshgrid(x0,x1)
z=theta0+theta1*x0+theta2*x1
#画3d图
ax.plot_surface(x0,x1,z)
ax.set_xlabel('area')
ax.set_ylabel('distance')
ax.set_zlabel("Monthly turnover")
plt.show()
结果：
2.与之前最小二乘法结果进行比较
最小二乘法：
多元线性回归方程为:y= 41.51X1+ -0.34 X2+ 65.32
梯度下降法：
多元线性回归方程为:y= 45.05X1+ -0.19 X2+ 5.37
发现结果不一致，但又相似。那我们就来看看原理！
3.梯度下降法原理
参考来源：https://www.jianshu.com/p/424b7b70df7b
一句话概括
梯度下降法（Gradient Descent，GD）是一种常用的求解无约束最优化问题的方法，在最优化、统计学以及机器学习等领域有着广泛的应用。
举例分析
假设这样一个场景：一个人需要从山的某处开始下山，尽快到达山底。在下山之前他需要确认两件事：
下山的方向
下山的距离
这是因为下山的路有很多，他必须利用一些信息，找到从该处开始最陡峭的方向下山，这样可以保证他尽快到达山底。此外，这座山最陡峭的方向并不是一成不变的，每当走过一段规定的距离，他必须停下来，重新利用现有信息找到新的最陡峭的方向。通过反复进行该过程，最终抵达山底。
在选择每次行动的距离时，如果所选择的距离过大，则有可能偏离最陡峭的方向，甚至已经到达了最低点却没有停下来，从而跨过最低点而不自知，一直无法到达山底；如果距离过小，则需要频繁寻找最陡峭的方向，会非常耗时。要知道，每次寻找最陡峭的方向是非常复杂的！同样的，梯度下降法也会面临这个问题，因此需要我们找到最佳的学习率，在不偏离方向的同时耗时最短。
梯度下降法的一般求解框架
1.给定待优化连续可微函数 J ( Θ ) J\left ( \Theta \right) J(Θ)、学习率 α \alpha α以及一组初始值 Θ 0 = ( θ 01 , θ 02 , ⋯ , θ 0 l , ) \Theta_{0}=\left ( \theta_{01}, \theta_{02}, \cdots, \theta_{0l}, \right ) Θ0=(θ01,θ02,⋯,θ0l,)
2.计算待优化函数梯度： ▽ J ( Θ 0 ) \triangledown J\left ( \Theta _{0} \right ) ▽J(Θ0)
3.更新迭代公式： Θ 0 + 1 = Θ 0 − α ▽ J ( Θ 0 ) \Theta^{0+1}=\Theta _{0}-\alpha \triangledown J\left ( \Theta _{0} \right ) Θ0+1=Θ0−α▽J(Θ0)
4.计算\Theta^{0+1}处函数梯度 ▽ J ( Θ 0 + 1 ) \triangledown J\left ( \Theta _{0+1} \right) ▽J(Θ0+1)
5.计算梯度向量的模来判断算法是否收敛： ∥ ▽ J ( Θ ) ∥ ⩽ ε \left\| \triangledown J\left ( \Theta \right ) \right \|\leqslant \varepsilon ∥▽J(Θ)∥⩽ε
6.若收敛，算法停止，否则根据迭代公式继续迭代
根据计算梯度时所用数据量不同，可以分为三种基本方法：批量梯度下降法、小批量梯度下降法以及随机梯度下降法。详情请见：https://www.jianshu.com/p/424b7b70df7b
简单实例
求： f ( x ) = f ( x 1 , x 2 ) = 1 3 x 1 2 + 1 2 x 2 2 f\left ( x \right )=f\left ( x_{1},x_{2}\right )=\frac{1}{3}x_{1}^{2}+\frac{1}{2}x_{2}^{2} f(x)=f(x1,x2)=31x12+21x22 函数的极小值点。
解：
设初始点为 x 1 = ( x 1 ( 1 ) , x 2 ( 1 ) ) ⊤ = ( 3 , 2 ) ⊤ x_{1}=\left ( x_{1}^{\left ( 1 \right )},x_{2}^{\left ( 1 \right )} \right )^{\top }=\left ( 3,2 \right )^{\top } x1=(x1(1),x2(1))⊤=(3,2)⊤，学习率设为 λ \lambda λ。
初始点处梯度为 g 1 = g ( x 1 ) = ( 2 , 2 ) ⊤ ≠ 0 g _{1}=g\left ( x _{1} \right )=\left ( 2,2 \right )^{\top }\neq 0 g1=g(x1)=(2,2)⊤=0，因此更新迭代公式带入原函数中，得：
f ( x 2 ) = f ( x 1 − λ g 1 ) = 10 3 λ 2 − 8 λ + 5 ， f\left( x_{2} \right )=f\left ( x_{1}-\lambda g_{1} \right )=\frac{10}{3}\lambda ^{2}-8\lambda +5， f(x2)=f(x1−λg1)=310λ2−8λ+5，此时 λ 1 ∗ = 6 5 \lambda _{1}^{*}=\frac{6}{5} λ1∗=56为函数极小点，因此：
x 2 = x 1 − λ 1 ∗ g 1 = ( 3 5 , − 2 5 ) ⊤ x_{2}=x_{1}-\lambda _{1}^{*}g_{1}=\left ( \frac{3}{5},-\frac{2}{5} \right )^{\top } x2=x1−λ1∗g1=(53,−52)⊤，一次迭代结束。
再将 x 2 x_{2} x2作为初始点，重复上面的迭代步骤，得到： x 3 = ( 3 5 2 , 2 5 2 ) ⊤ x_{3}=\left ( \frac{3}{5^{2}},\frac{2}{5^{2}} \right )^{\top } x3=(523,522)⊤。
根据规律显然可知，
x
k
=
(
3
5
k
−
1
,
(
−
1
)
k
−
1
2
5
k
−
1
)
⊤
x_{k}=\left(\frac{3}{5^{k-1}},\left ( -1 \right )^{k-1}\frac{2}{5^{k-1}} \right )^{\top }
xk=(5k−13,(−1)k−15k−12)⊤。
容易看出，本例中目标函数
f
(
x
)
f\left ( x \right )
f(x)是三维空间中的椭圆抛物面，其投影至二维空间上的等高线是一簇椭圆（如下图所示）。
f
(
x
)
f\left ( x \right )
f(x)的极小点就是这簇椭圆的中心
x
∗
=
(
0
,
0
)
⊤
x^{*}=\left ( 0,0 \right )^{\top }
x∗=(0,0)⊤。我们求得的迭代公式
{
x
k
}
\left \{ x_{k} \right \}
{xk}是逐渐趋近于x^{*}的，算法正确有效。
4.牛顿法原理
参考来源：https://www.jianshu.com/p/287b43e837c1
问题描述
f(x)是 R n R^n Rn上具有二阶连续偏导的函数，求 min x ∈ R n f ( x ) \min_{x\in R^n}f(x) minx∈Rnf(x)
原理
(1) 二阶泰勒展开
f
(
x
)
=
f
(
x
k
)
+
g
k
(
x
−
x
k
)
+
1
2
(
x
−
x
k
)
T
H
(
x
k
)
(
x
−
x
k
)
,
g
k
=
∇
f
(
x
)
f(x)=f(x_k)+g_k(x-x_k)+\frac12(x-x_k)^TH(x_k)(x-x_k),g_k=\nabla f(x)
f(x)=f(xk)+gk(x−xk)+21(x−xk)TH(xk)(x−xk),gk=∇f(x)
(2) 利用极值点导数为0求
x
k
+
1
x_{k+1}
xk+1
若
x
k
+
1
x_{k+1}
xk+1为极值点，则它的导数
g
k
+
1
=
0
g_{k+1}=0
gk+1=0，通过对二阶泰勒展开式求导可得：
g
k
+
1
=
g
k
+
H
k
(
x
k
+
1
−
x
k
)
=
0
⟹
x
k
+
1
=
x
k
−
H
−
1
g
k
g_{k+1}=g_k+H_k(x_{k+1}-x_k)=0\ \implies \ x_{k+1}=x_k-H^{-1}g_k
gk+1=gk+Hk(xk+1−xk)=0 ⟹ xk+1=xk−H−1gk
算法流程
输入:目标函数
f
(
x
)
f(x)
f(x)，梯度
g
(
x
)
=
∇
f
(
x
)
g(x)=\nabla f(x)
g(x)=∇f(x)，黑塞矩阵
H
(
x
)
H(x)
H(x)，精度要求
ϵ
\epsilon
ϵ
输出:
f
(
x
)
f(x)
f(x)的极小值点
(1) 取初值
x
0
x_0
x0，置
k
=
0
k=0
k=0
(2) 计算
g
k
g_k
gk
(3) 若
g
k
<
ϵ
g_k \lt \epsilon
gk<ϵ，则停止计算，得解
x
∗
=
x
k
x^*=x_k
x∗=xk
(4) 计算
H
k
=
H
(
x
k
)
,
求
p
k
,
H
k
p
k
=
−
g
k
H_k=H(x_k),求p_k,H_kp_k=-g_k
Hk=H(xk),求pk,Hkpk=−gk
(5) 计算
x
k
+
1
=
x
k
+
p
k
x_{k+1}=x_k+p_k
xk+1=xk+pk
(6)$ k=k+1$，转(2)
624
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
