梯度下降代码实现

import matplotlib.pyplot as plt
from sklearn.datasets.samples_generator  import make_blobs
import numpy as np
import argparse
def sigmoid_activation(x):
    return 1./(1+np.exp(-x))
ap = argparse.ArgumentParser()
ap.add_argument("-e","--epochs",type=float,default=100,help="# of epochs")
ap.add_argument("-a",'--alpha',type = float,default = 0.01,help = "learning rate")
args = vars(ap.parse_args())
(X,y) = make_blobs(n_samples = 250,n_features = 2,centers = 2,cluster_std = 1.05,random_state = 20)
# 表达式为1w0+w1x1+w2x2
X = np.c_[np.ones((X.shape[0])),X]
print('[INFO] starting training...')
# 初始化权重theta
theta = np.random.uniform(size=(X.shape[1],))
lossHistory = []
# 梯度下降训练
for epoch in np.arange(0,args["epochs"]):
    preds = sigmoid_activation(X.dot(theta))
    error = preds - y
    loss = np.sum(error**2)
    lossHistory.append(loss)
    print("[INFO epoch #{}],loss={:.7f}".format(epoch+1,loss))
    gradient = X.T.dot(error)/X.shape[0]
    theta += -args['alpha']*gradient
# 评估训练结果
for i in np.random.choice(250,10):
    activation = sigmoid_activation(X[i].dot(theta))
    label = 0 if activation < 0.5 else 1
    print('activation={:.4f}:predicted_label={},true_label={}'.format(activation,
                                                                     label,y[i]))
# 可视化结果
Y = (-theta[0]-(theta[1]*X))/theta[2]
fig = plt.figure()
plt.scatter(X[:,1],X[:,2],marker = 'o',c = y)
plt.plot(X,Y,'r')
plt.xlabel("Epoch #") 
plt.ylabel("Loss")
plt.show()

上述代码中激活函数的梯度如何求的?
对于二分类问题我们可以看做，在已知给定条件$$\mathbf{x}$$,假设已知变量$$\mathbf{w}$$时求预测$$\mathbf{y}$$最大的概率，$$h_{\theta}(x)=\frac{1}{1+e^{-\theta^Tx}}$$,因为sigmoid将函数映射到了0-1范围内，符合概率分布的特性，于是我们可以的出:

$$ \begin{equation} P(y=1|x;\theta) = h{\theta}(x)\ P(y=0|x;\theta) = 1-h{\theta}(x)\ \end{equation}

$$

$$P(y=1|x;\theta) = h_{\theta}(x)\$$表示在已知$$\mathbf{x}$$,参数$$\theta$$条件下判断$$\mathbf{x}$$为1的概率。将公式合起来:

$$ p(y|x;\theta) = (h{\theta}(x))^y(1-h{\theta}(x))^{1-y}

$$

假设样本$$\mathbf{x}$$之间相互独立，对于m组数据他们的联合概率密度函数为:

$$ \begin{equation} L(\theta)=p(\vec{y}|X;\theta)=\prod {i=1}^m(h{\theta})^{yi}(1-h{\theta}(x_i))^{(1-y_i)} \end{equation}

$$

我们的目标是找到一组$$\theta$$使得$$L(\theta)$$最小，即：求$$L(\theta)$$函数的最小值，对$$L(\theta)$$函数取自然对数不改变函数的最小值。于是得到:

$$ l(\theta)=\sum{i=1}^my_iln(h{\theta})+(1-yi)ln(1-h{\theta}(x_i))

$$

对这个函数求导数

$$ \begin{align} \frac{\partial l(\theta)}{\partial \thetaj} &= \frac{\partial l(\theta)}{\partial g(\theta^Tx)}\frac{\partial g(\theta^Tx)}{\partial \theta^Tx}\frac{\partial \theta^Tx}{\partial \theta_j}\ &=\left(\frac{y_i}{g(\theta^Tx)}-\frac{1-y_i}{1-h{\theta}(x_i)}\right) \cdot\left[g(\theta^Tx)(1-g(\theta^Tx))\right]\cdot x_i\ &=(y_i-g(\theta^x))x_i \end{align}

$$

对上面$$f(x) = \frac{1}{1+e^{-x}}$$求导可得:

$$ f'(x) = f(x)(1-f(x))

$$

所以上面的梯度等于error和x的矩阵乘法。

为什么求解的时候$$Y = (-theta[0]-(theta[1]*X))/theta[2]$$？j计算最佳拟合需要sigmod函数为0或者1，则$$\theta_0+\theta_1X_1+\theta_2X_2 =0$$则给定$$\theta_0$$的最佳拟合求出的$$X_2$$可以通过上面的方程求出。 源代码gradient_sk.py