机器学习之逻辑回归

一、逻辑回归的形式

如何使用线性模型进行分类任务？
可以用广义线性模型解决：求一个单调可微的函数将分类任务的真实标记与线性回归模型的预测值联系起来。

考虑二分类任务，其输出标记为{1，0}，而线性回归的预测值$z=w^{T}x+b$是实值，于是我们需要将实值转换为{0,1}值，最理想的是“单位阶跃函数”，但是它存在一个问题，那就是函数不连续，于是又找到一个替代函数sigmoid

如上图所示，这就是sigmoid的图形，当x>0时，y>0.5,x<0时，y<0.5。它的表达形式为：
$$y=\frac{1}{1+e^{-z}}$$

上式还可以变换为：
$$y=\frac{1}{1+e^{-(W^T+b)}}$$

$$W^T+b=ln\frac{y}{1-y}$$

若将y视为样本X作为正例的可能性，则1-y是其反例的可能性，两者的比值$\frac{y}{1-y}$称为“几率”，反映了x作为正例的可能性，对几率取对数则得到“对数几率”

$ln\frac{y}{1-y}$

由上式可看出，实际上是用线性回归模型的预测值去逼近真实标记的对数几率，因此，其对应模型称为“对数几率回归”，名字虽然是回归，但实际上做得是分类任务。
对数几率回归的优点：直接对分类可能性建模，无需事先假设数据分布，这样就避免了假设分布所带来的的不准确问题；它不是仅预测出“类别”，而是可得到近似概率预测，这对许多利用概率辅助决策的任务很有用。

二、参数的估计

如何确定w和b，若将y视为类后验概率估计$p(y=1|x)$,则上式可重写为

$ln\frac{p(y=1|x)}{p(y=0|x)}=w^{T}x+b$

$p(y=1|x)=\frac{1}{1+e^{-(w^{T}x+b)}}=\frac{e^{w^{T}x+b}}{1+e^{w^{T}x+b}}$

$p(y=0|x)=1-p(y=1|x)=\frac{1}{1+e^{w^{T}x+b}}$

其中 x ϵ R n , y ϵ { 0 , 1 } w ϵ R n x\epsilon R^n,y \epsilon \left \{0,1 \right \} w\epsilon R^n xϵRn,yϵ{0,1}wϵRn是权值向量，$bϵR$是偏置。

可将权值向量和特征向量加以扩充，即增广权值向量$\hat{W}=(w^1,w^2...,w^n,b{)}^T$,增广特征向量$\hat{x}=(x^1,x^2,...,x^n,1{)}^T$，则逻辑回归模型为：

$p(y=1|x)=\frac{e^{\hat{w}\hat{x}}}{1+e^{\hat{w}\hat{x}}}$

$p(y=0|x)=\frac{1}{1+e^{\hat{w}\hat{x}}}$

现在模型已经确定，就要估计参数了。

二分类对数几率回归参数学习——“最大似然估计”

给定训练数据集 D = { ( x 1 ^ , y 1 ) , ( x 2 ^ , y 2 ) , . . . ( x N ^ , y N ) } ， x i ^ = ( x 1 , x 2 , . . . , x n , 1 ) , y i ϵ { 0 , 1 } , i = 1 , 2 , 3 , . . . , n 。 D=\left \{ (\hat{x_1},y_1),(\hat{x_2},y_2),...(\hat{x_N},y_N)\right\}，\hat{x_i}=(x^1,x^2,...,x^n,1),y_i \epsilon \left\{0,1\right\},i=1,2,3,...,n。 D={(x1​^​,y1​),(x2​^​,y2​),...(xN​^​,yN​)}，xi​^​=(x1,x2,...,xn,1),yi​ϵ{0,1},i=1,2,3,...,n。

设$p(y=1|x)=g(\hat{w}\cdot \hat{x}),p(y=0|x)=1-g(\hat{w}\cdot \hat{x})$

似然函数：
$L(\hat{W})={\prod }_{i=1}^{N}P(y_i|{\hat{x}}_i)$

$L(\hat{W})={\prod }_{i=1}^N[g(\hat{w}\cdot {\hat{x}}_i){]}^{y_i}[1-g(\hat{w}\cdot {\hat{x}}_i){]}^{1-y_i}$

对数似然函数：
$l(\hat{W})=logL(\hat{W})$

$l(\hat{W})=log{\prod }_{i=1}^N[g(\hat{w}\cdot {\hat{x}}_i){]}^{y_i}[1-g(\hat{w}\cdot {\hat{x}}_i){]}^{1-y_i}$

$l(\hat{W})={\sum }_{i=1}^N[y_i\ast log(g(\hat{w}\cdot {\hat{x}}_i))+(1-y_i)log(1-g(\hat{w}\cdot {\hat{x}}_i))]$

最大似然估计，求解$W^{\ast }$
${\hat{W}}^{\ast }=\underset{\hat{W}}{\mathrm{argmax}}l(\hat{W})$

化为最小对数损失
${\hat{W}}^{\ast }=\underset{\hat{W}}{\mathrm{argmin}}-l(\hat{W})$

令${\hat{y}}_i=g(\hat{w}\cdot \hat{x})$,对数似然函数对$\hat{W}$求导得：

$\frac{\partial l(\hat{W})}{\partial \hat{W}}=-{\sum }_{i=1}^N[y_i\frac{{\hat{y}}_i(1-{\hat{y}}_i)}{{\hat{y}}_i}{x}_i-(1-y_i)\frac{{\hat{y}}_i(1-{\hat{y}}_i)}{1-{\hat{y}}_i}{x}_i]$

$\frac{\partial l(\hat{W})}{\partial \hat{W}}=-{\sum }_{i=1}^N{x}_i[y_i(1-{\hat{y}}_i)-{\hat{y}}_i(1-y_i)]$

$\frac{\partial l(\hat{W})}{\partial \hat{W}}=-{\sum }_{i=1}^N{x}_i(y_i-{\hat{y}}_i)$
注意：${{\hat{y}}_i}^{\prime }={\hat{y}}_i(1-{\hat{y}}_i)$，因为sigmoid函数的导数就是这样。

最后一步，我们采取梯度下降法对权值向量进行迭代更新

$${\hat{W}}_{t+1}={\hat{W}}_t+\alpha \sum _{i=1}^N{x}_i(y_i-{\hat{y}}_i^{wt})$$

其中$\alpha$表示学习率，${\hat{y}}_i^{wt}$表示模型在参数为${\hat{W}}_t$时的预测输出。

这样，我们的模型就出来了
$$y=\frac{1}{1+e^{-(wx+b)}}$$
在预测新的数据时，只需把数据放进这个模型，根据输出值比较设定的阈值来判断输出的类别。例如阈值设为0.5，如果y的输出值大于0.5，则判别类型为1，如果小于0.5，则判别类型为0.

三、代码实现

import numpy as np
import matplotlib.pyplot as plt
from sklearn import datasets
%matplotlib inline

#加载数据
def load_data():
	data = datasets.load_iris()
	X = data.data[:100,:2]
	y = data.target[:100]
	return X,y

#配置增广X矩阵
def data_matrix(X):
	data_mat = []
	for d in X:
		data_mat.append([1,*d])
	return np.array(data_mat)

#sigmoid函数
def sigmoid(z):
	return 1/(1+np.exp(-z))

#拟合数据
def fit(X,y,alpha=0.01,maxiters=200):
	m = y.size
	theta = np.zeros((len(X[0]),1),dtype=np.float32)
	for iter in range(maxiters):
		for i in range(m):
			h = sigmoid(np.dot(X[i],theta))
			theta += alpha * np.transpose([X[i]]) * (y[i]-h)
	return theta

#画出决策边界
X,y = load_data()
X = data_matrix(X)
theta = fit(X,y)
xx = np.arange(4,8)
yy = -(theta[0]+theta[1]*xx)/theta[2]
plt.scatter(X[:50,1],X[:50,2])
plt.scatter(X[50:,1],X[50:,2])
plt.plot(xx,yy)

结果如下图