逻辑回归（Logistic Regression） 是目前各行业最常用的分类方法，属于线性回归的拓展。
特点：
该算法联系了样本的特征和样本发生概率($\hat p = f(x)$),在运算上由于 概率值 本身是一个数值，因此该方法分类方法被称为回归方法。算法最终得到样例的预测概率值 $\hat p$ 用于 分类问题。所以逻辑回归既可以看作是回归算法，也可以看作是分类算法，通常作为分类算法使用，只可以解决二分类问题。
$$ \hat y = \left\{ \begin{aligned} 0,\hat p \le 0.5; \\ 1,\hat p \ge 0.5 \end{aligned} \right. $$
1、逻辑回归的基本原理
逻辑回归源于线性回归的拓展，具体思路是将线性回归的样本输出预测值到样本概率值的一个转换过程。在线性回归$\hat y = f(x) = \theta^T \cdot x_b$中，有预测的样本标记输出值 $\hat y \in (-\infty,\infty)$。为了使样本输出值 $\hat y $ 转换为概率值 $\hat p \in [0,1]$，可通过$Sigmoid$ 函数转换 样本值到 样本概率 $\hat p = \sigma(\theta^T\cdot x_b)$。
关于Sigmoid 函数
$$\sigma (t) = \frac {1}{1+e^{-t}}$$
Sigmoid 曲线具有如下性质:
- 函数曲线最左端无限趋近于0，函数曲线最右端无限趋近于1；值域 $\sigma (t) \in (0,1)$，因此很方便地用于表达概率；
- $ \left \{ \begin{array}. p > 0.5, t>0; \\ p<0.5,t<0 \end{array} \right.$
使用 $Sigmoid$ 函数转换样本值为样本概率 $\hat p = \frac {1}{1+e^{-(\theta^T\cdot x_b)}}$，进一步基于样本概率进行二分类任务 $\hat y = \left\{ \begin{array}. 1, \hat p \ge 0.5; \\0,\hat p \le 0.5. \end{array}\right.$
1.2 逻辑回归的优化目标
对于给定的样本集$X$和样本分类$y$，找到参数$\theta$, 使得用$Sigmoid$ 转换可以最大程度获得样本数据集$X$对应的分类输出$y$。
1.2.1 目标函数的定义
对于二分类问题，样本输出标记取值为离散变量 $y\in \{0,1\}$。使用 $Sigmoid $框架来获取样本特征到样本概率 ，再实现样本离散值标记的转换的时候。由于原始样本标记是一个离散值，而逻辑回归输出标记是一个分类概率值，意味着不能通过与类似线性回归等算法建立预测值与标记值的误差函数作为优化目标。
逻辑回归的损失函数定义建立在离散的样本标记上，定义如下 惩罚函数 ：
$cost = \left\{ \begin{array}. \small {如果 y=1 , \hat p 越小，cost越大} \\ \small {如果 y=0 , \hat p 越大，cost越大}\end{array}\right.$
在进行逻辑回归的时候，如果样本 真值为1，当前$\theta$下的预测概率值$\hat p $越小，则有惩罚越大；如果样本 真值为0，当前$\theta$下的预测概率值$\hat p $越大，则有惩罚越大；从而回归的目标为取得 最小惩罚下 对应的 参数$\theta$值。
1.2.2 惩罚函数的选取
$cost = \left\{ \begin{array}. -log(\hat p) \ \ \ \ \ \ \ \ \ \ \ if \ \ \ y = 1 \\ -log(1 - \hat p) \ \ \ \ if \ \ \ y = 0 \end{array}\right.$
- 当$y = 1$时，使用 $ cost(\hat p) = -log(\hat p)$ 作为惩罚函数，$\hat p \in [0,1]$ ，当 $\hat p \to 0$ 的时候有 惩罚$cost(\hat p) $ 趋于正无穷；只有$\hat p \to 1$ 时惩罚$cost(\hat p) $ 才会变小。
- 当$y = 0$时，使用 $ cost(\hat p) = -log(1- \hat p)$ 作为惩罚函数，$\hat p \in [0,1]$ ，当 $\hat p \to 1$ 的时候有 惩罚$cost(\hat p) $ 趋于正无穷；只有$\hat p \to 0$ 时惩罚$cost(\hat p) $ 才会变小。
进一步将分条件的 cost 函数联合为一个函数:
$$cost(\hat p) = \left\{ \begin{array}. -log(\hat p) \ \ \ \ \ \ \ \ \ \ \ if \ \ \ y = 1 \\ -log(1 - \hat p) \ \ \ \ if \ \ \ y = 0 \end{array}\right. \rightarrow cost(\hat p) = y*(-log(\hat p)) + (1-y)*(-log(1-\hat p ))$$
当$y =1 ,cost$ 函数仅表示为 $-log(\hat p) \ \ \ \ \ \ $的结果；
当$y =0 ,cost$ 函数仅表示为 $-log(1- \hat p)$的结果。全数据集X的样本损失函数表达为:
$$J(\theta) = \frac {1}{m} \sum _{i=1}^{m}{cost(\hat p)} = -\frac {1}{m} \sum _{i=1}^{m}{y^{(i)}log(\hat p^{(i)})+(1-y^{(i)}) log(1- \hat p^{(i)})}$$
$$\small {\bf 第(i) 个样本的概率估计值} \ \ \ \ \ \hat p^{(i)} = \sigma (\theta^Tx_b^{(i)}) = \frac {1}{1+e^{-\theta^Tx_b^{(i)}}}$$
回归目标是找到一组$\theta$是的$J(\theta)$得到最小值；逻辑回归的损失函数 无正规方程解，$J(\theta)$ 本质是一个 凸函数(Convex)，存在全局最优解，只能用 梯度下降法求解 。
1.2.3 $J(\theta)$的梯度
对于损失函数$J(\theta)$有如下形式：
$$J(\theta) = -\frac {1}{m} \sum _{i=1}^{m}{y^{(i)}log(\sigma(t))+(1-y^{(i)}) log(1- \sigma(t))} $$
- 根据求导链式法则对前半部分 $y^{(i)}log(\sigma(t))$ 求导
首先，令 $t = \theta^Tx_b^{(i)}$
$\because \sigma (t) = \frac {1}{1+e^{(-t)}} = (1+e^{(-t)})^{-1}$
$\therefore \sigma (t)' = -1(1+e^{(-t)})^{-2}*e^{(-t)}*(\ln{e^{-1}}) = (1+e^{-t})^{-2} \cdot e^{-t}$
$\because \log (\sigma {(t)})' = \frac {1}{\sigma (t)} \cdot \sigma {(t)}' = \frac {1}{\sigma (t)} \cdot (1+e^{-t})^{-2} \cdot e^{-t} = (1+e^{-t}) \cdot (1+e^{-t})^{-2} \cdot e^{-t} = (1+e^{-t})^{-1} \cdot e^{-t} = \frac {1+e^{-t} -1}{1+e^{-t}} = 1 - \frac {1}{1+e^{-t}}$
$\therefore \log (\sigma {(t)})' = 1- \sigma(t)$$\therefore[y^{(i)}log(\sigma(t))]' = \frac {\partial (y^{(i)} \log(\sigma(\theta^Tx_b^{(i)})))}{\partial(\theta_j)} = y^{(i)}*(1-\sigma(\theta^Tx_b^{(i)}))*x_{b\_j}^{(i)}$
..........................................................................................
- 根据求导链式法则对后半部分 $(1-y^{(i)})log(1 \sigma(t))$ 求导
$[\log(1-\sigma(t))]' = \frac {1}{1-\sigma(t)}*(-\sigma(t))' = \frac {1}{1-\sigma(t)}* -(1+e^{-t})^{-2} \cdot e^{-t} = -\frac {1+e^{-t}}{e^{-t}}*(1+e^{-t})^{-2} \cdot e^{-t} = -(1+e^{-t})^{-1} = -\sigma(t)$$\therefore [(1-y^{(i)})\log(1-\sigma(t))]' = \frac {\partial ((1-y^{(i)})\log(1-\sigma(\theta^Tx_b^{(i)})))}{\partial \theta_j} = (1-y^{(i)})*(-\sigma(\theta^Tx_b^{(i)}))*x_{b\_j}^{(i)}$
..........................................................................................
联合前后部分函数的求导结果为 $\nabla J(\theta_j)$
$\nabla J(\theta_j) = -\frac {1}{m} \sum _{i=1}^{m} y^{(i)}*(1-\sigma(\theta^Tx_b^{(i)}))*x_{b\_j}^{(i)} + (1-y^{(i)})*(-\sigma(\theta^Tx_b^{(i)}))*x_{b\_j}^{(i)}$
$\nabla J(\theta_j) = \frac {1}{m} \sum _{i=1}^{m} (\sigma(\theta^Tx_b^{(i)}) - y^{(i)} )\cdot x_{b\_j}^{(i)} = \frac {1}{m} \sum _{i=1}^{m} (\hat y^{(i)} - y^{(i)} )\cdot x_{b\_j}^{(i)}$
推广到 所有$\theta$
$\nabla J(\theta)= \begin{bmatrix} \partial J/ \partial \theta _{0} \\ \partial J/ \partial \theta _{1} \\ ... \\ \partial J/ \partial \theta _{n} \end{bmatrix} = \frac {1}{m} \begin{bmatrix} \sum_{i=1}^{m}{\sigma ( X_{b}^{(i)}\theta - y^{i}) \cdot ( \ 1 \ )} \\ \sum_{i=1}^{m}{\sigma ( X_{b}^{(i)}\theta - y^{i}) \cdot (x^{i}_1)} \\... \\ \sum_{i=1}^{m}{\sigma ( X_{b}^{(i)}\theta - y^{i}) \cdot (x^{i}_n)} \end{bmatrix} = \frac {1}{m} \cdot X^T_b\cdot (\sigma(X_b\theta) - y)$逻辑回归损失函数的梯度结果与线性回归极其类似！！！
1.2.4 参数$\theta$ 梯度下降法求解一般过程
### Prepare data
import numpy as np
from sklearn import datasets
iris = datasets.load_iris()
from sklearn.model_selection import train_test_split
Train_X,Test_X,Train_Y,Test_Y = train_test_split(iris.data[iris.target < 2,:2],iris.target[iris.target < 2],test_size= 0.2,random_state=666)
X_b = np.hstack([np.ones((len(Train_X),1)),Train_X])
y = Train_Y
### sigmoid function
def _sigmoid(t):
return 1 / (np.exp(-t) + 1)
### Loss Function
def J(theta,X_b,y):
y_hat = _sigmoid(X_b.dot(theta))
try:
return - np.sum(y*np.log(y_hat)+(1-y)*np.log(1-y_hat))/len(y)
except:
return float('inf') ### 异常检测，如果出现异常返回浮点数最大值float('inf')
### Derivative function
def dJ(X_b,y,theta):
return X_b.T.dot(_sigmoid(X_b.dot(theta)) - y) / len(X_b)
### Gradient Descent
eta = .000001 ### 学习率
epsilon = 1e-8 ### 定义一个误差范围，当初始点开始滚动，如果当前点位与一点位函数值误差小于定义的epsilon误差范围，认为当前点为就是最优解。
theta_lst = [np.random.randint(10,size = X_b.shape[1])+ np.random.normal(size = X_b.shape[1])] ### 随机化滚落初始点
while True:
gradient = dJ(X_b,y,theta_lst[-1])
theta_lst.append(theta_lst[-1] - eta * gradient)
if(abs(J(theta_lst[-1],X_b,y) - J(theta_lst[-2],X_b,y)) < epsilon):break
if(len(theta_lst) > 1e4): break ### 循环1e4次以后如果没有取得最优解将退出循环
print(theta_lst[-1])
1.2.5 scikit-learn 框架下线性分布数据的逻辑回归
import numpy as np
from sklearn import datasets
iris = datasets.load_iris()
from sklearn.model_selection import train_test_split
Train_X,Test_X,Train_Y,Test_Y = train_test_split(iris.data[iris.target < 2,:2],iris.target[iris.target < 2],random_state=666)
from sklearn.linear_model import LogisticRegression
log_reg = LogisticRegression()
log_reg.fit(Train_X,Train_Y)
log_reg.score(Test_X,Test_Y)
log_reg.decision_function(Test_X) ### 获取 θ^T·x_b 的决策值 ，
2、非线性数据的逻辑回归
对于 线性分布的数据，意味着其样本在特征空间的类别分布可以被 一条直线或(超平面) 所分割，对应 决策边界 表示为线性方程 $\theta^Tx_b = 0$。对于 非线性分布的数据 ， 其在特征空间中的类别分布无法被 一条直线或(超平面) 所分割。
在 多项式回归 一章中，非线性数据的回归问题可以通过构造特征转换为线性回归问题。同理在 逻辑回归 中处理非线性分布数据的判别问题中也可以通过 添加多项式项（构造特征） 再进行 逻辑回归，从而 获得任意形状的决策边界。
2.2 scikit-learn 框架下非线性分布数据的逻辑回归
### Prepare data
import numpy as np
import matplotlib.pyplot as plt
np.random.seed(666)
X = np.random.normal(0,1,size = (200,2))
y = np.array(X[:,0]**2 +X[:,1] < 1.5 ,dtype='int') ### 令曲线 -x_0^2 - x_1 + 1.5 = 0 下方的样本分类为 1
for _ in range(20): ### 生成分类噪音
y[np.random.randint(200)] = 1
### Darw sample distribution
plt.scatter(X[y == 0 ,0],X[y == 0 ,1])
plt.scatter(X[y == 1 ,0],X[y == 1 ,1])
plt.show()
### Logistic Regression
from sklearn.preprocessing import PolynomialFeatures
from sklearn.pipeline import Pipeline
from sklearn.preprocessing import StandardScaler
from sklearn.linear_model import LogisticRegression
def PolynomialLogisticRegression(degree):
return Pipeline([
('poly',PolynomialFeatures(degree=degree)),
('std_scaler',StandardScaler()),
('logistic_reg',LogisticRegression())
])
poly_log_reg = PolynomialLogisticRegression(degree= 2)
poly_log_reg.fit(Train_X,Train_Y)
poly_log_reg.score(Test_X,Test_Y)
3、逻辑回归中添加正则化
对于线性算法，当添加上多项式项时模型复杂度随着随着多项式项阶数(degree)的上升而指数级增大，进而容易发生对 训练集过拟合 的问题，在逻辑回归中也是如此(表现在决策边界的规整情况)。解决模型过拟合问题的常规处理为向模型添加 正则化 ，通过约束参数大小从而提高模型的泛化能力降低方差。
3.2 逻辑回归的正则项
- 损失函数添加 L1 正则项
$$C\cdot J(\theta) + L_1 = C\cdot J(\theta) + \sum_{i=1}^{n}{|\theta_i|}$$ - 向损失函数添加 L2 正则项
$$C\cdot J(\theta) + \alpha L_2 = C\cdot J(\theta) + \frac {1}{2}\sum_{i=1}^{n}{\theta_i}^2$$
3.3 scikit-learn 框架下的逻辑回归添加正则项
from sklearn.preprocessing import PolynomialFeatures
from sklearn.pipeline import Pipeline
from sklearn.preprocessing import StandardScaler
from sklearn.linear_model import LogisticRegression
def PolynomialLogisticRegression(degree,c,penalty = 'l2'):
return Pipeline([
('poly',PolynomialFeatures(degree=degree)),
('std_scaler',StandardScaler()),
('logistic_reg',LogisticRegression(solver= 'liblinear',C = c,penalty = penalty)) ### 参数C 指定正则化约束强度; 参数penalty 指定添加正则项的类别(L1,L2)
])
poly_log_reg = PolynomialLogisticRegression(degree= 20,c = .1,penalty='l1',)
poly_log_reg.fit(Train_X,Train_Y)
poly_log_reg.score(Test_X,Test_Y)
4、逻辑回归对多分类问题的解决方案
4.1 OvR (One vs Rest)
将多个类别每次选取一个类(One)设定为 1，其余类别((Rest))则设定为 0 ，从而转换为 多个二分类任务。这样对于 $n$ 个类别就会创建 $n$ 个分类模型，对新样例进行 $n$ 个模型的预测，选择分类得分最高的作为新样例的分类。
缺点：时间复杂度高，效率低。
sklearn的逻辑回归 OvR 多分类支持
from sklearn.linear_model import LogisticRegression
OvR_logreg = LogisticRegression(multi_class='ovr')
4.2 OvO (One vs One)
将多个类别每次挑选其中两个类别构建一个二分类模型，这样对于数据集中的 $n$ 个类别就会组合出$C_{n}^2$ 个 二分类模型，然后分别应用这些分类器对新样例进行预测，最后选择这 $C_{n}^2$ 个分类器对 新样例 的预测结果中出现最多的类别作为新样例的分类。
优点：由于每个分类器都是对两个真实类比的比较，不像 OvR 任务中 Rest 部分存在混杂类别带来的干扰信息，因此 OvO 分类器更倾向于得到新样例的真实类别，也即分类准确度更高。
缺点：时间复杂度比 OvR 还高，效率低。
sklearn的逻辑回归 OvO 多分类支持
from sklearn.linear_model import LogisticRegression
OvO_logreg = LogisticRegression(multi_class='multinomial',solver="newton-cg")
4.3 sklearn 框架下转换二分类机为多分类器
### 创建逻辑回归二分类器
from sklearn.linear_model import LogisticRegression
log_reg = LogisticRegression()
from sklearn.multiclass import OneVsOneClassifier
ovo = OneVsOneClassifier(log_reg) ### 封装逻辑回归二分类器为 ovo 多分类器
ovo.fit((Train_X,Train_Y)
ovo.score(Test_X,Test_Y)
from sklearn.multiclass import OneVsRestClassifier
ovr = OneVsRestClassifier(log_reg) ### 封装逻辑回归二分类器为 ovr 多分类器
ovr.fit((Train_X,Train_Y)
ovr.score(Test_X,Test_Y)
