逻辑回归
1 Logistic回归分析基本思想
Logistic回归实际上是一种分类方法，称为概率型非线性回归模型，是研究二分类观察结果与其影响因素之间关系的一种多变量分析方法。二分类问题，其输出结果为\(y\)属于{0,1}，而线性回归模型产生的预测值\(z=wx+b\)是实值，所以我们要把实值\(z\)转化为0/1值。最理想的是“单位阶跃函数”。 “回归”的意思就是要找到最佳拟合参数，其中涉及的数学原理和步骤如下：
1. 找到一个合适的分类函数实现分类，比如单位阶跃函数、Sigmoid函数。
2. 构建损失函数cost来测量预测值与实际值的偏差，为使得回归拟合最佳，就要满足偏差最小。
3. 记\(J(ω)\)表示回归系数为\(ω\)时的偏差，那么求最佳回归参数\(ω\)就转换成了求\(J(ω)\)的最大值。采用梯度上升法求解。
1.1 分类函数
我们想要的函数应该是，能接受所有的输入然后预测出类别。例如，在两个类的情况下上述函数输出0或1。Sigmoid函数可以满足这个要求，Sigmoid函数具体的计算公式如下：
\[\sigma(z)=\frac{1}{1+e^{-z}}\]
下图给出了Sigmoid函数在不同坐标尺度下的两条曲线图。当x为0时，Sigmoid函数值为0.5。随着x的增大，对应的Sigmoid值将逼近于1；而随着x的减小，Sigmoid值将逼近于0。
z <- c(-100:100)
y <- 1/(1+exp(-z))
plot(z,y,type = 'l',main="逻辑回归曲线",col = "red",pch = 22)
abline(h=0.5,v=0)
因此，为了实现Logistic回归分类器，我们可以在每个特征上都乘以一个回归系数，然后把所有的结果值相加为\(Z\)，即
\[Z=g(x)=ω_0x_0+ω_1x_1+ω_2x_2+…+ω_nx_n=ω^TX\]
向量\(x\)是特征变量，是输入数据，向量\(w\)是回归系数。将这个总和\(Z\)作为Sigmoid函数的输入值，进而得到一个范围在0~1之间的数值。任何大于0.5的数据被分入1类，小于0.5即被归入0类。所以，Logistic回归也可以被看成是一种概率估计。
接下来详细谈下logistic回归模型:
设条件概率\(P(y=1|x)=p\)为根据观测量相对于某事件发生的概率，那么logistic回归模型可以表示为: \[P(y=1|x)=\pi(x)=\frac{1}{1+e^{-g(x)}}\] 那么在\(x\)条件下\(y\)不发生的概率为:
\[P(y=0|x)=1-\pi(x)=\frac{1}{1+e^{g(x)}}\] 所以事件发生与不发生的概率之比为: \[\frac{P(y=1|x)}{P(y=0|x)}=\frac{p}{1-p}=e^{g(x)}\] 这个比值称为事件的发生比，简记为odds。对odds取对数得到: \[ln\frac{p}{1-p}=g(x)=ω_0x_0+ω_1x_1+ω_2x_2+…+ω_nx_n=ω^TX\]
可以看出Logistic回归都是围绕一个Logistic函数来展开的。
那接下来的问题是确定函数z中最佳的回归系数w。
1.2 cost函数
那么怎么找到最佳的回归系数，这里是构造cost函数来最小化拟合回归方程的误差。
对于特征变量的线性形式： \[Z=g(x)=ω_0x_0+ω_1x_1+ω_2x_2+…+ω_nx_n=ω^TX\] 分类函数： \[h(z)=\frac{1}{1+e^{-z}}\]
因此进一步推导得： \[h(z)=h_w(x)=\frac{1}{1+e^{-w^Tx}}\] 而由1.1推导得，对于任意确定的x和w，有： \[P(y=1|x,w)=h_w(x)\\P(y=0|x,w)=1-h_w(x)\] 当\(x>0\) 时，y=1,即：（可结合逻辑回归图理解这里） \[1/2<h_w(x)<1\\P(y=1|x,w)=h_w(x)>1/2\] 当\(x<0\) 时，y=0,即：（可结合逻辑回归图理解这里） \[0<h_w(x)<1/2\\P(y=0|x,w)=1-h_w(x)>1/2\] 所以，无论如何取值，\(P(y|x,w)\)都大于等于1/2,该概率越大，越接近1，表示落入某分类概率越大，那么分类越准确，预测值与实际值差异就越小。 所以\(P(y|x,w)\)可以表示预测值与实际值的差异，且\(P(y|x,w)\)越大表示差异越小。
将\(P(y|x,w)\)归纳为： \[P(y|x,w)=h_w(x)^y(1-h_w(x))^{1-y}\] 对于每个点\((x^i,y^i)\)都有上述等式成立，那么为使得其成立的概率最大。这里用极大似然估计来计算参数。
取对数似然函数： \[J(w)=\sum^m_{i=1}y^ilogh_w(x^i)+(1-y^i)log(1-h_w(x^i))\]
所以似然函数\(J(w)\)越大，预测越准确，偏差越小。因此，\(J(w)\)可以来表示预测值与实际值的偏差，所以这里的\(J(w)\)可以看成cost函数，不过注意这个cost函数是当函数值越大时，其偏差越小。
2 Logistic回归参数估计——梯度下降法
上面已经得到了cost函数，但是我们还是不知道怎么去求解分类函数里面的线性分类器\(Z\)，即求解\(g(x)\)中的\(w\)向量。
这里介绍梯度下降算法来求解。
梯度下降法（gradient descent），又名最速下降法（steepest descent）是求解无约束最优化问题最常用的方法，它是一种迭代方法，每一步主要的操作是求解目标函数的梯度向量，将当前位置的负梯度方向作为搜索方向。梯度下降不一定能够找到全局的最优解，有可能是一个局部最优解。当然，如果损失函数是凸函数，梯度下降法得到的解就一定是全局最优解。
梯度下降法特点：越接近目标值，步长越小，下降速度越慢。
下面将通过公式来说明梯度下降法。
需要拟合求解的函数\(g(w)\) ,因为\(x\)是数据中已知的特征值: \[g(w)=ω_0x_0+ω_1x_1+ω_2x_2+…+ω_nx_n=ω^TX\] 接下来的目标是将该函数通过样本的拟合出来，得到最佳的函数模型。因此构建损失函数\(J(w)\)（目的是通过求解min \(J(w)\)，得到在最优解下的\(w\)向量），其中的每一项\(g(w,x^i)-y^i\)表示在已有的训练集上我们的拟合函数与y之间的残差，计算其平方损失函数作为我们构建的风险函数（这里采用最小二乘法构造损失函数，在逻辑回归中也可采用最大似然估计构造损失函数从而估计参数）。 \[min J(w)=\frac{1}{2}\sum^m_{i=0}(g(w,x^i)-y^i)^2\] 要使得最小\(J(w)\)，则\(J(w)\)对\(w\)求导等于零。 \[\nabla{J(w)}=\sum^m_{i=0}(g(w,x^i)-y^i)x^i\] 对于多个特征变量，这里的\(x^i\)为特征变量值矩阵，\(x^i_j\)为第\(j\)个特征的值的向量。
梯度下降法的目标：要找到\(\nabla{J(w)}\)等于零时\(w\)的值。
梯度下降法步骤：
1. 给定终止误差\(\sigma>0\)（比如设置为\(1×10^{-5}\)），令\(k:=0，w:=0\)；
2. 计算\(\nabla{J(w)}\), 若\(||\nabla{J(w)}||<\sigma\),停止迭代，输出\(w\),否则进行第3步；
3. 构造负梯度方向\(-\nabla{J(w)}\).迭代\(w\)向量，对每个\(w_j:=w_j-t\sum^m_{i=0}(g(w,x^i)-y^i)x^i_j\)，t为步长，一般取1，得到\(w\)向量后重新代入第2步，直到满足第2步中的条件后退出循环。得到\(w\)向量。
得到\(w\)向量后我们就求得了逻辑回归中的分类函数。当然在R中我们可以直接调用包来完成逻辑回归模型的构建。
例：以AER包中的Affairs数据集作为例子。该数据集是关于婚姻出轨，其中affairs变量表示出轨次数，数据集中还包括结婚时间、教育、宗教等其它变量。由于affairs为正整数，为了进行Logistic回归先要将其转化为二元变量。
data(Affairs,package = 'AER')
summary(Affairs)## affairs gender age yearsmarried children
## Min. : 0.000 female:315 Min. :17.50 Min. : 0.125 no :171
## 1st Qu.: 0.000 male :286 1st Qu.:27.00 1st Qu.: 4.000 yes:430
## Median : 0.000 Median :32.00 Median : 7.000
## Mean : 1.456 Mean :32.49 Mean : 8.178
## 3rd Qu.: 0.000 3rd Qu.:37.00 3rd Qu.:15.000
## Max. :12.000 Max. :57.00 Max. :15.000
## religiousness education occupation rating
## Min. :1.000 Min. : 9.00 Min. :1.000 Min. :1.000
## 1st Qu.:2.000 1st Qu.:14.00 1st Qu.:3.000 1st Qu.:3.000
## Median :3.000 Median :16.00 Median :5.000 Median :4.000
## Mean :3.116 Mean :16.17 Mean :4.195 Mean :3.932
## 3rd Qu.:4.000 3rd Qu.:18.00 3rd Qu.:6.000 3rd Qu.:5.000
## Max. :5.000 Max. :20.00 Max. :7.000 Max. :5.000Affairs$naffairs[Affairs$affairs==0]=0
Affairs$naffairs[Affairs$affairs>0]=1
Affairs$naffairs <- factor(Affairs$naffairs)
table(Affairs$naffairs)##
## 0 1
## 451 150可以看出naffairs这一新列中，出轨的例子有150起，未出轨的例子有451起。
lmodel <- glm(formula=naffairs~gender+age+yearsmarried+children+religiousness+education+occupation+rating,family=binomial(link='logit'),data=Affairs)
summary(lmodel)##
## Call:
## glm(formula = naffairs ~ gender + age + yearsmarried + children +
## religiousness + education + occupation + rating, family = binomial(link = "logit"),
## data = Affairs)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -1.5713 -0.7499 -0.5690 -0.2539 2.5191
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 1.37726 0.88776 1.551 0.120807
## gendermale 0.28029 0.23909 1.172 0.241083
## age -0.04426 0.01825 -2.425 0.015301 *
## yearsmarried 0.09477 0.03221 2.942 0.003262 **
## childrenyes 0.39767 0.29151 1.364 0.172508
## religiousness -0.32472 0.08975 -3.618 0.000297 ***
## education 0.02105 0.05051 0.417 0.676851
## occupation 0.03092 0.07178 0.431 0.666630
## rating -0.46845 0.09091 -5.153 2.56e-07 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 675.38 on 600 degrees of freedom
## Residual deviance: 609.51 on 592 degrees of freedom
## AIC: 627.51
##
## Number of Fisher Scoring iterations: 4通过多元逻辑回归结果看出是否出轨较大程度取决于age，yearsmarried，religiousness和rating这几个变量。
设age，yearsmarried，religiousness和rating这几个变量分别为\(x_1,x_2,x_3,x_4\)。
那么该分类函数形式为\[y=\frac{1}{1+e^{-(-0.04426x_1+0.09477x_2-0.32472x_3-0.46845x_4)}}\]
这个模型到底效果怎么样，还需要进行模型评估。这个留给同学们动手试试。
3 本章汇总
| 名称 | 类别 | 功能 |
|---|---|---|
| data() | 函数 | 加载R内置数据集 |
| summary() | 函数 | 打印模型拟合的结果 |
| table() | 函数 | 将数据以表格形式呈现 |
| factor() | 函数 | 转换为因子变量 |
| glm() | 函数 | 广义线性回归模型 |