机器学习算法之： 逻辑回归 logistic regression (LR)

by joey周琦

LR介绍

逻辑回归属于probabilistic discriminative model这一类的分类算法

probabilistic discriminative mode这类算法的思路如下：
- 直接建模P(Ck|x)
- 利用最大似然估计和训练数据，估计出模型中的参数

该类想法相对于生成模型（probabilistic generated model) 有参数较少的优点。因为生成模型需要P(x|Ck)和先验概率P(Ck).

LR是工业界最长用的分类算法之一，其主要原因，个人认为有几点如下：

• 训练速度快，扛得住大数据
• 模型可解释度、可理解程度高，根据每个特征的系数，就可以判断出该特征在模型中的重要性，帮助判断模型是否合理
• 可以接受的精度

本文，对LR做一个简单的总结

LR二分类

首先做下简单的符号说明，在下述推导中，N为样本个数，M为特征数目，x∈RM为特征向量， K为可分类的总数，P(Ck|x)表示在特征向量的前提下，判别为第k类的概率。w∈RM为LR模型的参数，即训练LR的主要目的就是要得到w。

对于K=2类，建模如下:

P(C1|x)=y(x)=σ(wTx)=11+e−wTx

σ是logistic sigmoid function, 它有个性质（1）：

dσ(a)da=σ(1−σ)

其中yn=σ(an)且 an=wTx。假设我们有数据集{ xn,tn }, n=1...N, N为数据集的大小，n表示第n个样本，tn=0or1表示第n个样本的真实类别，设t=(t1,…,tN),那么似然函数可以写为如下形式:

p(t|w)=∏n=1Nytnn(1−yn)1−tn

最大化似然函数等价于最大化对数似然函数可以写为

lnp(t|w)=∑n=1Ntnlnyn+(1−tn)ln(1−yn)

我们要最大化对数似然函数，就是等价于对数似然函数的相反数，则最小化目标函数可以写为：

J=−lnp(t|w)=−∑n=1Ntnlnyn+(1−tn)ln(1−yn)

目标函数对参数w求导，利用上述σ函数的性质1，可以推到得到：

∂J∂w=∑n=1N(yn−tn)xn

可以利用梯度下降法对参数w进行更新.更新公式为：

w=w−λ∂J∂w=w−λ∑n=1N(yn−tn)xn

其中λ是学习率，可以自己设置一个比较小的数字，或者通过其他算法计算（比如牛顿法）。可以看到每次更新都要有一个N想得加和，这样很需要时间，所以更新策略可以采用“随机梯度下降法”，每次只用一个样本，如下：

w=w−λ∂J∂w=w−λ(yn−tn)x

这样大大节约了训练时间，也不会特别影响精度。有些系统为了避免过拟合，目标函数中可以加入L1或者 L2正则项（这里采用L2)，可以避免||w||2过大,也就是为了避免模型为了拟合训练数据而变得过于复杂。加入正则项后，目标函数变为

J′=J+α2||w||2

其中α为正则项惩罚系数。那么随机梯度的更新策略也就变为了

w=w−λ∂J′∂w=w−λ(yn−tn)xn−λαw

LR多分类

关于多分类的，在特征x∈RM下，模型判别为k类的概率如下：

P(Ck|x)=yk(x)=exp(wTkx)∑Kj=1exp(wTjx)

，
其中wj∈RM表示第j个模型的系数.另外，tn∈RK表示的时第n个样本的真实分类, tn向量只有一个元素为1，其余为0，若tnk=1则，表示该样本实际是第k类。T=t1,t2,…,tN表示所有样本的真实分类

那么最大似然函数则可以表示为：

P(T|w1,…,wK)=∏n=1N∏k=1KP(Ck|xn)tkn

取负并取对数，则得到目标函数：

J=−logP(T|w1,…,wK)=−∑Nn=1∑Kk=1tnklogynk

可以推导（有点麻烦）：

∂J∂wj=∑Nn=1(ynj−tnj)xn

所以，梯度下降法对参数wj进行更新公式：

wj=wj−λ∑Nn=1(ynj−tnj)xn

随机梯度下降法：

wj=wj−λ(ynj−tnj)xn

随机梯度下降法+L2正则，目标函数变为：

J′=J+α2||w||2

则参数的更新变为：

wj=wj−λ(ynj−tnj)xn−αλwj

值得注意的是：

• 数据来自kaggle练习比赛的https://www.kaggle.com/c/digit-recognizer/data, 在这里我只用了train.csv数据，其中10%作为测试数据
• 由于数据的feature是像素，0~255，比较大，所以学习率要设置的低一些，或者将feature归一化。
• 加入正则项并不一定可以提高精度
• 结果如果2分类，区分数字0与9，正确率99%。如果是多分类，正确率为84%左右。

LR二分类

import numpy as npp_test = [ ]p_train = [ ]lines = [line.rstrip('\n') for line in open('train.csv')]label = 0features = np.zeros(784) #28*28 featurefor line in lines[1:]:    line = line[0:-1]    data = line.split(",")    for i in range(len(data)):        data[i]=float(data[i])    data=np.array(data)    if np.random.rand()<0.1:        p_test.append(data)    else:        p_train.append(data)#Train LR for 2 class learing_r = 0.001w = np.random.rand(784) - 0.5for img in p_train:    label = img[0]    value = np.array(img[1:])    if label==9: # classify 9 and 0 , ignore others        t=1    elif label==0:        t=0    else:        continue    y = 1.0/(1+np.exp(-np.dot(w,value)))    w = w - learing_r * (y-t) * value#Test LRright = 0 wrong = 0for img in p_test:    label = img[0]    value = np.array(img[1:])    if label==9: # classify 9 and 0 , ignore others        t=1    elif label==0:        t=0    else:        continue    y = 1.0/(1+np.exp(-np.dot(w,value)))    if y>0.5:        pt=1    else:        pt=0    if pt==t:        right+=1    else:        wrong+=1print right*1.0/(right+wrong)

LR多分类

'''Created on 2015 8 23@author: joeyqzhou'''import numpy as npimport bigfloatp_test = [ ]p_train = [ ]lines = [line.rstrip('\n') for line in open('train.csv')]label = 0features = np.zeros(784) #28*28 featurefor line in lines[1:]:    line = line[0:-1]    data = line.split(",")    for i in range(len(data)):        data[i]=float(data[i])    data=np.array(data)    if np.random.rand()<0.1:        p_test.append(data)    else:        p_train.append(data)#Train LR for multiple class learing_r = 0.000001W = [ ]K = 10 # 10 class for i in range(K):    W.append(np.random.rand(784)*0.1 - 0.05)for img in p_train:    label = img[0]    value = np.array(img[1:])    y = np.zeros(K)    temp_y = np.zeros(K) #unnormalized y    for i in range(K):        temp_y[i] = bigfloat.exp( np.dot(W[i], value) )    for i in range(K):        y[i] = temp_y[i]/np.sum(temp_y)    for i in range(K): #update the parameter        if abs(i-label) < 0.0001:            W[i] = W[i] - learing_r * (y[i]-1) * value        else:            W[i] = W[i] - learing_r * (y[i]-0) * value#Test LR for multiple classright_count = 0whole_count = 0for img in p_test:    label = img[0]    value = np.array(img[1:])    y = np.zeros(K)    for i in range(K):        y[i] =  np.dot(W[i], value)  #not necessary to normalize     inx = np.argmax(y)    if abs(inx - label)<0.001:        right_count += 1    whole_count += 1print right_count*1.0/whole_count