线性回归

• 回归分析
• 目标函数：线性回归方程 $y = wx + b$
• 一个或多个自变量和因变量之间的关系进行建模(其中$\theta_i$为权重，$\theta_0$为bias偏置值)：
  • 一维特征：$h_\theta(x)=\theta_0+\theta_1x_1$
  • 二维特征：$h_\theta(x)=\theta_0+\theta_1x_1+\theta_2x_2$
  • N维特征：$\sum_{i=0}^N\theta_ix_i=\theta^TX$

梯度下降算法

构建损失函数

• 未得到目标线性方程，得到上述公式中的$\theta^T$参数
• 为确定所选的$\theta_T$效果好坏使用损失函数（loss function）来评估$h(x)$函数的好坏
• 损失函数如下(MSE误差平方和）：
  $J(\theta)=\frac{1}{2m}\sum_{i=1}^m(h(\theta^Tx_i)-y_i)^2 \tag{1}$

梯度下降法

• 使用了导数的概念，对点$x_0$的导数反映了函数在该点处的瞬时变化速率
  $f^\prime(x_0)=\lim_{\vartriangle x \to \infty} =\frac{\vartriangle y}{\vartriangle x}=\lim_{\vartriangle x \to \infty} =\frac{f(x_0+\vartriangle x)-f(x_0)}{\vartriangle x} \tag{2}$
• 推广到多维函数中，梯度反映了多维图形中变化速率最快的方向
• 起始时导数为正，$\theta$减小后并以新的$\theta$为基点重新求导，一直迭代就会找到最小的$J(\theta)$
  
• StochasticGradientDescent(SGD)随机梯度下降
  $\text{Repeat} \{ \\ \dotsb\\ \theta_j = \theta_j -\alpha((h(\theta^Tx_i)-y_i)x_i)\\ \dotsb\\ \}$
  其中：$x_i$为随机的样本点
• 代码实现

import random
import matplotlib.pyplot as pltx = [(0.,3), (1.,3) ,(2.,3), (3.,2), (4.,4), (0.,3) , (1.,3.1) ,(2.,3.5), (3.,2.1) , (4.,4.2)]
#y[i] is the output of y = theta0 * x[0] + theta1 * x[1] +theta2 * x[2]
y = [95.364,97.217205,75.195834,60.105519,49.342380, 100.364,100.217205,100.195834,100.105519,12.342380]epsilon = 0.001
#learning rate
alpha = 0.001
diff = [0,0]
error1 = 0
error0 =0
m = len(x)#init the parameters to zero
theta0 = 0
theta1 = 0
theta2 = 0
epoch = 0error_array = []
epoch_array = []
while True:#calculate the parameters# 线性回归：h(x) = theta0  + theta1 * x[i][0] + theta2 * x[i][1] # 损失函数(评测指标MSE)：累和 (1/2) *  (y - h(x)) ^ 2# d((1/2) *  (y - h(x)) ^ 2) / dtheta0 = (y - h(x)) * (-1) * (1) # theta0 = theta0 - (-alpha * (y - h(x))* 1 )# theta1 = theta1 - (-alpha * (y - h(x))* x[i][0])# theta2 = theta2 - (-alpha * (y - h(x))* x[i][1])# 1. 随机梯度下降算法在迭代的时候，每迭代一个新的样本，就会更新一次所有的theta参数。i = random.randint(0, m - 1)# (y - h(x))diff[0] = y[i]-( theta0 * 1 + theta1 * x[i][0] + theta2 * x[i][1] )# - (y - h(x))xgradient0 = - diff[0]* 1gradient1 = - diff[0]* x[i][0]gradient2 = - diff[0]* x[i][1]# theta = theta - (  - alpha * (y - h(x))x )theta0 = theta0 - alpha * gradient0theta1 = theta1 - alpha * gradient1theta2 = theta2 - alpha * gradient2#theta3#calculate the cost functionerror1 = 0# 此处error为一个相对的Error值。for i in range(len(x)):error1 += (y[i]-( theta0 * 1 + theta1 * x[i][0] + theta2 * x[i][1]))**2 error1 = error1 / m    print("delta  error {}".format(abs(error1-error0)))error_array.append(error1)epoch_array.append(epoch)if epoch == 500:breakelse:error0 = error1epoch += 1print(' theta0 : %f, theta1 : %f, theta2 : %f, sgd error1 : %f, epoch : % f'%(theta0,theta1,theta2,error1, epoch))print('Done: theta0 : %f, theta1 : %f, theta2 : %f'%(theta0,theta1,theta2))plt.plot(epoch_array, error_array, color='blue', linewidth=3)
plt.xlabel('epochs')
plt.ylabel('errors')
plt.show()

• BatchGradientDescent(BGD)批量随机梯度下降
  $\text{Repeat} \{ \\ \dotsb\\ \theta_j = \theta_j -\alpha(\frac{1}{m}\sum_{i=0}^m(h(\theta^Tx_i)-y_i)x_i)\\ \dotsb\\ \}$
• 代码实现

import matplotlib.pyplot as plt#Training data set
#each element in x represents (x0,x1,x2)
x = [(0.,3) , (1.,3) ,(2.,3), (3.,2) , (4.,4), (0.,3) , (1.,3.1) ,(2.,3.5), (3.,2.1) , (4.,4.2)]
#y[i] is the output of y = theta0 * x[0] + theta1 * x[1] +theta2 * x[2]
y = [95.364,97.217205,75.195834,60.105519,49.342380, 100.364,100.217205,100.195834,100.105519,12.342380]#learning rate
alpha = 0.001 
diff = [0,0]
error1 = 0
error0 =0
m = len(x)#init the parameters to zero
theta0 = 0
theta1 = 0
theta2 = 0
sum0 = 0
sum1 = 0
sum2 = 0epoch = 0
error_array = []
epoch_array = []
while True:#calculate the parameters# 线性回归：hi(x) = theta0 + theta1 * x[i][1] + theta2 * x[i][2]  # 损失函数：(1/2) 累加 * (y - h(x)) ^ 2# theta = theta - 累和(  - alpha * (y - h(x))x )# 1. 随机梯度下降算法在迭代的时候，每迭代一个新的样本，就会更新一次所有的theta参数。#calculate the parameters# 2. 批梯度下降算法在迭代的时候，是完成所有样本的迭代后才会去更新一次theta参数print(m)for i in range(m):#begin batch gradient descentdiff[0] = y[i]-( theta0 + theta1 * x[i][0] + theta2 * x[i][1] )sum0 = sum0 + ( -alpha * diff[0]* 1)sum1 = sum1 + ( -alpha * diff[0]* x[i][0])sum2 = sum2 + ( -alpha * diff[0]* x[i][1])#end  batch gradient descenttheta0 = theta0 - sum0 / mtheta1 = theta1 - sum1 / mtheta2 = theta2 - sum2 / msum0 = 0sum1 = 0sum2 = 0#calculate the cost functionerror1 = 0for i in range(m):error1 += ( y[i]-( theta0 + theta1 * x[i][0] + theta2 * x[i][1] ) )**2error1 = error1 / merror_array.append(error1)epoch_array.append(epoch)if epoch == 200:breakelse:error0 = error1epoch += 1print(' theta0 : %f, theta1 : %f, theta2 : %f, bgd error1 : %f, epoch: %f'%(theta0,theta1,theta2,error1,epoch))print('Done: theta0 : %f, theta1 : %f, theta2 : %f'%(theta0,theta1,theta2))
plt.plot(epoch_array, error_array, color='blue', linewidth=3)
plt.xlabel('epochs')
plt.ylabel('errors')
plt.show()

对比两类梯度下降法各有利弊，SGD的计算量小，但随机性强，有时可以快速收敛，有时则收敛的很慢，loss曲线抖动下降；MGD的计算量大，但包含批量的样本，收敛效果比较好，loss曲线平滑下降。

Logistic Regression算法

• 针对分类问题，利用回归问题的思路，解决分类问题。

sigmoid函数

• sigmoid函数是一个s形的曲线
• 取值在(0, 1)之间
• 在远离0的地方函数的值会很快接近0或1
• $g(z)=\frac{1}{1+e^{-z}}$

构造目标函数

• $\theta_0+\theta_1x_1+\dotsc+\theta_nx_n=\sum_{i=0}^n\theta_ix_i=\theta^TX$
• $h_\theta(x)=g(\theta^TX)=\frac{1}{1+e^{-\theta^TX}}$
• $p(y=1|x,\theta)=h_\theta(x) \\p(y=0|x,\theta)=1-h_\theta(x) \tag{3}$
• 将预测结果转化成预测某类别的概率公式(3)中y为真实的类别

构造损失函数-极大似然估计

• 根据上述公式构建损失函数，目的：通过$\theta^T$的取值使预测正确类别的概率最大化
• $p(y|x,\theta)=(h_\theta(x))^y(1-h_\theta(x))^{1-y}$
• Lossfunction: $L(\theta)=\prod_{i=0}^mp(y_i|x_i,\theta)=\prod_{i=0}^m(h_\theta(x_i))^{y_i}(1-h_\theta(x_i))^{1-y_i}$
• 由于累乘的式子难以求导，因此对上式进行对数处理，转化为累加
• $l(\theta)=\log L(\theta)=\sum_{i=0}^m(y_i(\log h_\theta(x_i)+(1-y_i)\log (1-h_\theta)))$
• 转化为类似线性回归问题中求最小损失函数的问题：
  • $J(\theta)=-\frac{1}{m}l(\theta)$

梯度下降

• 首先求解损失函数梯度
  $\begin{aligned}\frac{\partial}{\partial \theta_j}J(\theta) & = -\frac{1}{m}\sum_{i=1}^m\left[y_i\frac{1}{h_\theta(x_i)}\frac{\partial}{\partial \theta_j}h_\theta(x_i)-(1-y_i)\frac{1}{1-h_\theta(x_i)}\frac{\partial}{\partial \theta_j}h_\theta(x_i)\right] \\& = \frac{1}{m} \sum_{i=1}^m\left(h_\theta(x_i)-y_i\right)x_i^j\\\end{aligned}$
• 使用梯度下降求解最优参数
  $\theta_j = \theta_j - \alpha\frac{1}{m} \sum_{i=1}^m\left(h_\theta(x_i)-y_i\right)x_i^j$

多分类问题

• 解决方案：one vs rest
  • 适合场景：类别可以同时出现，互相不互斥的情况
  • 参考资料：one-vs-rest与one-vs-one以及sklearn的实现

优化算法：牛顿法

• 牛顿法比梯度下降法收敛的要快（迭代次数更少）
  • 梯度下降法每次只从当前位置选一个上升速度最大的方向走一步
  • 牛顿法在选择方向时，不仅会考虑上升速度是否够大，还会考虑你走了一步之后，上升速度是否会变得更大（即二阶导数，类似物理上的加速度）
• 红色的牛顿法的迭代路径，绿色的是梯度下降法的迭代路径。
  
• 由于最速梯度下降收敛速度并不“最速”，局部搜索中最速下降的方向和全局的最小值方向并不一致，所以就有后来改进的方法，包括牛顿法以及拟牛顿法。

切线法

• 牛顿法的几何意义
  
• 从上图可以总结出如下公式：
  $\theta_j = \theta_j - \frac{J'(\theta_j)}{J''(\theta_j)}$
  其中分母上的二阶导函数可以理解为正则项，当一阶导数下降方向不是最优方向，进行修改，而且由于二阶导数包含数值，因此也不用设置学习率，下降的步长是最优的。

另一种理解方式

• 随时函数最小值 = 损失函数导数为0
• 将损失函数进行泰勒展开：
  $f(x) \approx f\left(x^{(k)}\right) + g^T_k\left(x-x^{(k)}\right)+\frac{1}{2}\left(x-x^{(k)}\right)^TH\left(x^{(k)}\right)\left(x-x^{(k)}\right) \tag{4}$
  其中$g(x)$为一阶导数，$H$为二阶导数
• 对(4)式等式两端进行求导，转化为梯度下降：
  $\triangledown f(x)=g_k +H_k \left( x-x^{(k)}\right)\\ g_k +H_k \left( x-x^{(k)}\right) = 0 \\ x^{(x+1)} = x^{(x)} - H_k^{-1}g_k$

  最终的结果与切线法得到的结果表达式是一样的

改进：拟牛顿法

• 牛顿法中需要计算$H_k$二阶导数矩阵，其尺寸为$N\times N$的，其中N为特征的数量，计算代价很大
• 而拟牛顿法不在使用$H_k$二阶导数矩阵，而是选取相似矩阵

Softmax Regression算法

• Softmax解决多分类的分类器，即类标签y的取值大于等于2
• 类标记为: $y(i) \in \left\{ 1,2,\cdots,k \right\}$
• 假设函数为对于每一个样本估计其所属的类别的概率$p(y=j∣x)$
• 每一个样本估计所属类别概率为：
  $p\left(y^{(i)}|x^{(x)};\theta\right)=\frac{e^{\theta^T_jx^{(i)}}}{\sum_{l=1}^ke^{\theta^T_lx^{(i)}}}$
  即对预测数的各类别的值取$log$后除以所有类别概率值取$log$的总和

Softmax回归代价函数

• 类似于Logistic回归，在Softmax的代价函数中引入指示函数$I\{⋅\}$，其具体形式为：
  $I\left \{ \text{expression}\right\} = \begin{cases}0 \quad \text{if}\ \text{expression}=\text{false}\\1 \quad \text{if}\ \text{expression}=\text{true} \end{cases}$
• 那么，对于Softmax回归的代价函数为交叉熵：
  $J(\theta)=-\frac{1}{m}\left[\sum_{i=1}^m\sum_{j=1}^kI\{y^{i}=j \}\log \frac{e^{\theta^T_jx^{(i)}}}{\sum_{l=1}^ke^{\theta^T_lx^{(i)}}} \right]$
• Sofrmax回归求解
  • 对上述的代价函数，可以使用梯度下降法进行求解，首先对其进行求梯度：
    $\begin{aligned}\triangledown_{\theta_j}J(\theta) & =-\frac{1}{m}\sum_{i=1}^m\left[\triangledown_{\theta_j}\sum_{j=1}^kI\left\{y^{(i)}=j \right\}\log \frac{e^{\theta^T_jx^{(i)}}}{\sum_{l=1}^ke^{\theta^T_lx^{(i)}}} \right]\\& =-\frac{1}{m}\sum_{i=1}^m\left[I\left\{y^{(i)}=j \right\} \frac{\sum_{l=1}^ke^{\theta^T_lx^{(i)}}}{e^{\theta^T_jx^{(i)}}}\cdot \frac{e^{\theta^T_jx^{(i)}}\cdot x^{(i)}\cdot \sum_{l=1}^ke^{\theta^T_lx^{(i)}}-e^{\theta^T_jx^{(i)}}\cdot x^{(i)}\cdot e^{\theta^T_jx^{(i)}}}{\left( \sum_{l=1}^ke^{\theta^T_lx^{(i)}}\right)^2}\right]\\ & = -\frac{1}{m}\sum_{i=1}^m \left[ I\left\{y^{(i)}=j \right\} \cdot \frac{\sum_{l=1}^ke^{\theta^T_lx^{(i)}}- e^{\theta^T_jx^{(i)}}}{\sum_{l=1}^ke^{\theta^T_lx^{(i)}}} \cdot x^{(i)}\right] \end{aligned}$
    $\theta_j = \theta_j - \triangledown_{\theta_j}J(\theta)$

L1/L2正则化

• 为了防止过拟合，拟合过程中通常都倾向于让权值尽可能小，最后构造一个所有参数都比较小的模型
• 参数值小的模型比较简单，能适应不同的数据集，也在一定程度上避免了过拟合现象
• 对于一个线性回归方程，若参数很大，那么只要数据偏移一点点，就会对结果
  造成很大的影响
• 但如果参数足够小，数据偏移得多一点也不会对结果造成什么影响

L1

• L1正则化是指权值向量$w$中各个元素的绝对值之和，通常表示为$\left\|w \right\|_1$
• 如下图所示，图中等值线是J0的等值线，黑色方形是L函数的图形。在图中，当J0等值线与L首次相交的地方就是最优解。上图中J0与L在L的一个顶点处相交，这个顶点就是最优解。
  
• 正则化和“带约束的目标函数”是等价的对L1正则化建立数学模型如下:
  $J(\theta)=\sum_{i=1}^n\left(y_i - \theta^T x_i \right)^2\\ s.t.\quad \left\|\theta \right\|_1\leqslant m$
• 通过拉格朗日乘子法，可以变为如下形式：
  $J(\theta)=\sum_{i=1}^n\left(y_i - \theta^T x_i \right)^2 + \lambda \left(\left\|\theta \right\|_1-m\right)$

L2

• L2正则化是指权值向量$w$中各个元素的平方和然后再求平方根（在Ridge回归的L2正则化项有平方符号），通常表示为$\left\|w \right\|_2$
• 如下图所示，二维平面下L2正则化的函数图形是个圆，与方形相比，被磨去了棱角。
  
• 对比L1和L2两图，采用L1范数时平方误差项等值线与正则化项等值线的交点出现在坐标轴上，即$w_1$或$w_2$为0，而采用L2范数时，两者的交点常出现在某个象限中，即$w_1$或$w_2$均非0；换而言之，采用L1范数比L2范数更容易于得到稀疏解。
• 正则化和“带约束的目标函数”是等价的对L1正则化建立数学模型如下:
  $J(\theta)=\sum_{i=1}^n\left(y_i - \theta^T x_i \right)^2\\ s.t.\quad \left\|\theta \right\|_2^2\leqslant m$
• 通过拉格朗日乘子法，可以变为如下形式：
  $J(\theta)=\sum_{i=1}^n\left(y_i - \theta^T x_i \right)^2 + \lambda \left(\left\|\theta \right\|_2^2-m\right)$

L1和L2对比

其中c为超参数，控制L1和L2的限制区域大小，c越小对损失函数的约束越大；图中像素点表示[0, 1]的取值，0为白色1为黑色。可以看出L1下是很多权值为0，做到了稀疏性。

参考资料：L1与L2范数

正则化目的

• 预防或修复过拟合情况
• 其他方式：
  • 数据层面：
    • 增加数据（行上扩展）
    • 减少feature，特征选择
  • 模型层面：
    • 模型融合

Ridge与Lasso

• L1正则化和L2正则化可以看做是损失函数的惩罚项
• Lasso回归
  • 对于线性回归模型，使用L1正则化的模型建叫做Lasso回归
  • Lasso回归的损失函数，式中加号后面一项即为L1正则化项
• Ridge
  • Ridge回归的损失函数，式中加号后面一项$\alpha \left\|\theta \right\|_2^2$即为L2正则化项

ElasticNet

对于ElasticNet回归，可以看做是Lasso和ridge回归的组合
$J(\theta)=\sum_{i=1}^n\left( y_i - \theta^Tx_i\right)^2+\lambda_1(\left\|\theta \right\|_1)+\lambda_2(\left\|\theta \right\|_2^2)\\ \min_\theta J(\theta)$