1 线性目标的梯度优化
损失函数：
算法1 : 批量梯度下降BGD
每次迭代使用所有样本来对参数进行更新。
损失函数：
代数形式：
矩阵形式：
更新：
代数形式伪代码：
矩阵形式伪代码：
算法2：随机梯度下降SGD
每次迭代使用一个样本来对参数进行更新。
一个样本的损失函数：
代数形式伪代码：
矩阵形式伪代码
算法3：小批量梯度下降法
每次迭代使用n_batch个样本对参数进行更新
代数形式伪代码：
矩阵形式伪代码
1. 线性目标的sgd优化，写出伪代码。mini_batch怎么更新？
2. SGD并行实现？
=============================
# 梯度优化代码
import numpy as np import pandas as pdfrom sklearn.model_selection import train_test_split from sklearn.linear_model import LinearRegression #基准模型 data_df = pd.read_csv('kc_house_data.csv') train_data, test_data = train_test_split( data_df, train_size=0.8, random_state=1234 ) # sklearn模型 X_train = train_data['sqft_living'].values.reshape(-1,1) y_train = train_data['price'].values.reshape(-1,1) X_test = test_data['sqft_living'].values.reshape(-1,1) y_test = test_data['price'].values.reshape(-1,1) model = LinearRegression(normalize=True).fit(X_train, y_train ) y_pred = model.predict(X_test) print('intercept:',model.intercept_) print('coel:', model.coef_) print('MSE: %.2f' % np.mean((y_pred - y_test) ** 2)) # 模型 def linear(X,w): X = np.hstack((np.ones((X.shape[0],1)), X)) y = np.dot(X,w) return y # 损失函数 def cost_function( X, y, w ): m = X.shape[0] e = np.dot(X,w) - y J = 0.5/m * np.dot(e.T, e) return J[0] # 标准化 def normalize( X ): mu = np.mean(X,axis =0) sigma = np.std(X,axis = 0) X = (X-mu)/sigma return X,mu,sigma # 均方误差 def MSE( X, y, w ): return np.mean((linear(X,w) - y) ** 2) #----------------------------- # 批量梯度下降法 def batch_gradient_descent( X, y, eta=0.1, eps=1e-4, num_iters=2000 ): converged = False J = np.zeros( num_iters ) m = X.shape[0] i_iters = 0 w = np.zeros((X.shape[1],1)) while not converged and i_iters < num_iters: grad_w = np.dot( X.T, X.dot(w) - y ) w -= eta*grad_w/m dist_grad = np.sqrt(np.dot(grad_w.T,grad_w)) J[i_iters] = cost_function( X, y, w ) print('Iteration:', i_iters+1, ': J(w):', J[i_iters], "||J(w)||:", dist_grad) if dist_grad < eps: cur= True i_iters += 1 return w X, mu, sigma = normalize( X_train ) X = np.hstack( (np.ones((X_train.shape[0],1)), X) ) w = batch_gradient_descent( X, y, eta = 0.1, eps=1e-4 ) print ('Intercept:', w[0] - np.dot(mu/sigma,w[1:])) print( 'Slope :', w[1:]/sigma[:,np.newaxis]) print ("MSE : %.2f" % MSE( (X_test-mu)/sigma, y_test, w )) #---------------------------------------------- # 随机梯度下降法 def stochastic_gradient_descent( X, y, eta=0.1, eps=1e-4, num_iters=2000 ): cur = False J = np.zeros( num_iters ) n_batch = 1 ## i_iters = 0 w = np.zeros((X.shape[1],1)) while not cur and i_iters < num_iters: batch = np.random.choice(X.shape[0], n_batch) X_batch, y_batch = X[batch], y[batch] grad_w = np.dot( X_batch.T, X_batch.dot(w) - y_batch ) w -= eta*grad_w/n_batch # dist_grad = np.sqrt(np.dot(grad_w.T,grad_w)) J[i_iters] = cost_function( X_batch, y_batch, w ) print('Iteration:', i_iters+1, ': J(w):', J[i_iters], "||J(w)||:", dist_grad) if dist_grad < eps: cur = True i_iters += 1 return w X, mu, sigma = normalize( X_train ) X = np.hstack( (np.ones((X_train.shape[0],1)), X) ) w = stochastic_gradient_descent( X, y, eta = 0.1, eps=1e-4 ) print ('Intercept:', w[0] - np.dot(mu/sigma,w[1:])) print( 'Slope :', w[1:]/sigma[:,np.newaxis]) print ("MSE : %.2f" % MSE( (X_test-mu)/sigma, y_test, w )) # --------------------------------------- # 小批量随机梯度下降法 def mini_gradient_descent( X, y, eta=0.1, eps=1e-4, num_iters=2000 ): cur = False J = np.zeros( num_iters ) n_batch = 5000 ## i_iters = 0 w = np.zeros((X.shape[1],1)) while not cur and i_iters < num_iters: batch = np.random.choice(X.shape[0], n_batch) X_batch, y_batch = X[batch], y[batch] grad_w = np.dot( X_batch.T, X_batch.dot(w) - y_batch ) w -= eta*grad_w/n_batch # n_batch为所选择的批量个数 dist_grad = np.sqrt(np.dot(grad_w.T,grad_w)) J[i_iters] = cost_function( X_batch, y_batch, w ) print('Iteration:', i_iters+1, ': J(w):', J[i_iters], "||J(w)||:", dist_grad) if dist_grad < eps: cur = True i_iters += 1 return w X, mu, sigma = normalize( X_train ) X = np.hstack( (np.ones((X_train.shape[0],1)), X) ) w = mini_gradient_descent( X, y, eta = 0.1, eps=1e-4 ) print ('Intercept:', w[0] - np.dot(mu/sigma,w[1:])) print( 'Slope :', w[1:]/sigma[:,np.newaxis]) print ("MSE : %.2f" % MSE( (X_test-mu)/sigma, y_test, w ))
代码写的再多，不如自己尝试一遍，很多东西都是自己操作了才学得比较快。这里的数据以及大部分知识点多是从公众号学习，这里就不具体说明。很多东西都是要自己找到的才会真正舍得去学。数据文件在网盘，http://pan.baidu.com/share/link?shareid=1375610667&uk=209896182。可以自己尝试下载数据练习，里面也有各种jupyter，还有岭回归和拉索回归的，感兴趣可以自行去学习。
2 LR回归原理推导
假设数据服从伯努利分布，通过极大似然函数的方法，用梯度下降来求解参数，来达到将数据分类，一般是二分类。
2.1 对率回归推导
1. sigmoid函数推导
逻辑回归解决的是二分类问题，首先通过训练模型来学习x和y的映射关系，有了映射关系就容易对x和y进行预测分类。如何定义这种映射关系，这里假定y的取值为0或者1(也可以是其他类)，两种可能，自然而然服从伯努利分布，这种映射关系就可以通过条件概率p(y|x)来表示。p(y=1|x)和p(y=0|x)，对于二分类就可设定一个阈值，0.5。
如何通过条件概率来描述x和y之间的关系呢?对率回归是线性模型，但是无法直接表示成
，那么就可以通过广义线性模型来实现。一般广义线性模型的一般形式， 
，其中g(y)为单调函数，g(y)就可被作为联系函数。（为什么选择sigmoid函数，一般先考虑的是单位阶跃函数将任意的实数转换为0/1的概率值，用它来当成联系函数来判断属于哪个类。但是单位阶跃函数的一个缺点就是在零点不连续也不单调，而联系函数需要单调连续。）而sigmoid能够将(-无穷，+无穷)的值域映射到(0,1)，这样就可以得到合理的概率值，而且单调连续，可以作为联系函数。
这样就得到了
，将
代入，得到
用条件概率表示就是
得到
2. 广义线性模型推导
指数族分布：是一类分布的总称，它的概率密度函数一般形式是：
其中，
称为该分布的自然参数；T(y)为充分统计量，通常为y本身；
为配分函数，保证概率表达式加和为1，保证式子是合格的概率密度函数；b(y)是关于随机变量y的函数。常见伯努利和正态均为指数族分布。
证明：伯努利分布是指数族分布？
化成指数族的一般形式：
对应指数族分布的一般形式，
广义线性模型三假设：
1. 给定x的条件下，假设随机变量y服从指数族分布
2. 给定x条件下，目标是得到一个模型h(x)能预测出T(y)的期望值。
3. 假设该指数族分布中的自然参数和x呈线性关系，即
满足这三条假设即为广义线性模型。
对数几率回归是在对一个二分类进行建模，并假设被建模的随机变量y取值为0或1，可以很自然地假设y服从伯努利分布。如果想要构建一个线性模型来预测在给定x的条件下y取值的话，可以考虑使用广义线性模型来建模。
假设1：伯努利分布服从指数族分布。
假设2：
，得到
假设3：
2.2 损失函数推导
损失函数： 
由式子，
假设数据都是独立的，由伯努利分布，得到一个样本发生的条件概率表示为：
进而得到它们一起发生的概率（似然函数）：
取对数，取负，得到损失函数
为什么用似然函数？
我们的目标是预测为某一类的出现概率最大，每个样本预测都要得到最大的概率。前面我们得到了一个样本的条件概率，极大似然估计就会将所有的样本考虑进去，来使得观测到的样本的出现概率最大，所以有了累乘形式的似然函数。 （为什么累乘？）基于条件独立的假设，总的条件概率就可以表示为每一个条件概率的累乘形式。
为什么用极大似然函数(交叉熵)作为损失函数而不用最小二乘（欧氏距离，均方损失）？
均方损失：假设误差是正态分布，适合线性的输出（回归问题），特点是对于与真实值差别越大，惩罚力度越大，不适合分类问题。
交叉熵损失：假设误差是伯努利分布，可以视为预测概率分布与真实概率分布的相似程度。多应用在分类的问题。
均方误差对参数的偏导的结果都乘以了 sigmoid的导数，而sigmoid导数在其变量值很大或很小时趋近于0，所以均方误差的偏导可能接近于0。
而参数更新公式： 参数 = 参数 – 学习率 * 损失函数对参数的偏导
在偏导很小时，参数更新速度会变得很慢，而在接近于0时，参数几乎不更新，出现梯度消失的情况。反观交叉熵对参数的偏导就没有sigmoid导数，所以不存在梯度消失的问题。
均方损失：
使用梯度下降法的条件损失函数时凸函数。而对最小二乘的损失函数求导，
可以知道，J(w)对w不是凸函数，不能用代价函数。
交叉熵：
当越接近于1，越接近于0。预测值与真实值完全相同，其损失函数为0。
当越接近于0，越接近于。
当越接近于1，越接近于。
当越接近于0，越接近于0。预测值与真实值完全相同，其损失函数为0。
写成矩阵形式：
2.3 梯度下降法
代数形式伪代码：
矩阵形式伪代码：
==================================
# LR代码
import numpy as np #随机生成样本数据，二分类问题，每类样本生成5000个 np.random.seed(43) num = 5000 x1 = np.random.multivariate_normal([0,0],[[1,0.75],[0.75,1]],num) #d多元正态分布矩阵，[0,0]为均值，然后是正定对称半正定矩阵。 x2 = np.random.multivariate_normal([1, 4], [[1, .75],[.75, 1]], num) X = np.vstack((x1, x2)).astype(np.float32) X = np.hstack( (np.ones((X.shape[0],1)), X) ) y = np.hstack((np.zeros(num),np.ones(num))) # 模型 def sigmoid(x): return 1 / (1 + np.exp(-x)) # 损失函数 def log_likelihood(X, y, w): # 首先按照标签来提取正样本和负样本的下标 pos, neg = np.where(y==1), np.where(y==0) pos_sum = np.sum(np.log(sigmoid(np.dot(X[pos], w)))) neg_sum = np.sum(np.log(1-sigmoid(np.dot(X[neg], w)))) return -(pos_sum + neg_sum) # 返回cross entropy loss # 实现逻辑回归模型 def logistic_regression(X, y, num_steps, eta): m = X.shape[1] w = np.zeros(X.shape[1]) for step in range(num_steps): error = sigmoid(np.dot(X,w)) - y grad_w = np.matmul(X.T, error) w = w - eta/m * grad_w # 每隔一段时间，计算一下损失函数，看看有没有变化 # 正常情况下， 它会慢慢变小，最后收敛 if step % 10000 == 0: print (log_likelihood(X, y, w)) return w w = logistic_regression(X, y, num_steps = 100000, eta = 5e-5) print ("参数w, b分别为: ", w[1:],w[0])# sklearn from sklearn.linear_model import LogisticRegression # C设置一个很大的值，意味着不想加入正则项 clf = LogisticRegression(fit_intercept=True, C = 1e15) clf.fit(X, y) print ("(sklearn)逻辑回归的参数w, b分别为: ", clf.coef_, clf.intercept_, )
=============================================================================================
写在后，这是第一次写博客，平常喜欢在word或txt记录，也习惯了word的公式，就不打latex了，公式直接截图了。发出来的初衷有两个，一个是怕万一本地的文档像上个电脑突然去世，里面的资料都丢失了。还一个也是希望，对上述公式或者知识点是以个人的见解来总结或者直接拿过来，希望能够看到大家的看法，以及纠正，在这也对直接拿过来的这些原出处没办法做引用说明说声抱歉。
也是在自学的路上，更新会比较慢，这也是对自己的一种鞭挞吧。
后续，还会有想对牛顿法，拟牛顿法进行补充，在能够认为掌握之后。
愿正在阅读的你，以及仍热爱的你，共勉。
