《机器学习实战》logistic回归：关于’此处略去了一个简单的数学推导‘的个人理解

正在看《机器学习实战》这本书的朋友，在看到logistic回归的地方，可能会对P78页的梯度上升算法代码以及P79这里的这句话弄的一头雾水：“此处略去了一个简单的数学推导，我把它留给有兴趣的读者”。这句话就是针对下面这段我贴出来的代码中的gradAscent函数说的。

我当时是想不明白为什么特征值矩阵的转置：dataMatrix.transpose()，它与真实类别/预测类别二者之间的差做矩阵乘法得到的结果，就被直接当作梯度了。

于是上网搜了下关于这部分内容的解释，baidu搜到的前两条链接里的内容都是基于均方误差（最小二乘法）作为代价函数对于此处做的数学解释。其他的搜索结果也未给予详细的数学推导过程。

而看过Andrew.N.G公开课logistic回归这一章的朋友，想必不会忘记吴老师关于logistic回归的代价函数处的额外补充，有别于线性回归，logistic回归如果依旧采用均方误差作代价函数，那么应用梯度下降法进行优化效果不会很理想（由于非凸，存在局部最优值），再翻看一下李航老师的《统计学习方法》一书，也能发现logistic回归的“策略”也是求取对数似然函数的极大值为目标，而非求取均方误差函数的极小值。

二来，上述代码是“梯度上升算法”的实现，如果拿均方误差函数作为代价函数，那岂不是助纣为虐，代价越来越大么？

基于上述两点疑问，决定自己做一下数学推导，推导过程如下：


推导过程除了符号复杂一些，真的不难。。。就是单纯的求偏导数而已

接下来是对代码中的过程进行下数学验证：


到这里，我们就可以发现，正常的数学推导得到的梯度跟代码中dataMatrix.transpose()*error是一模一样的。

补充：作者也是刚刚入门机器学习，本着较真的劲儿尝试着用均方误差做代价函数，看一看是不是存在很多局部最优值。。。
以这本书附带的数据集做为例子，设置了10000个权重向量，观察一下代价函数的图像。
以下是代码：

# -*- coding: utf-8 -*-
"""
Created on Wed Sep 26 19:40:20 2018

@author: _Miracle
本文件为验证logistic回归使用均方误差函数作代价函数不可行
"""

import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

def loadData():
    dataMat = []
    labelMat = []
    fr = open(r'C:\Users\_Miracle\Desktop\machinelearninginaction-master\Ch05\testSet.txt')
    for line in fr.readlines():
        lineArr = line.strip().split()
        dataMat.append([float(lineArr[0]),float(lineArr[1])])
        labelMat.append(int(lineArr[2]))
    return dataMat,labelMat

def sigmoid(intX):
    return 1.0/(1+np.exp(-intX))


sampleData,classLabel = loadData()
sampleDataMat = np.array(sampleData)

print ('这是特征数据值(100X2),100个样本，每个样本两个特征')
print (sampleDataMat)

#比较笨的方法，生成二维普通的list，然后调用
m,n = np.shape(sampleDataMat)
newListofList = []
for i in range(m):
    newList = []
    newList.append(classLabel[i])
    newListofList.append(newList)
classLabelMat = np.array(newListofList)

print ('这是每个样本的类别(100X1)，100个样本就有100个实际类别')
print (classLabelMat)

a = np.arange(m)
a = np.tile(a,(1,m))
#将a变为形如[[0，1，2，0，1，2，···，0，1，2]]重复arange(100)100次的ndarray
b = np.arange(m)
b = np.tile(b,(m,1))
b = b.transpose()
b = b.reshape(1,10000)
#将b处理成形如[[0,0,0,1,1,1,2,2,2]]这种0-99各重复100次的ndarray
#将a，b处理成这个样子是为了生成一个2x10000的矩阵，使得有10000种权重搭配的可能
weights = np.vstack((a,b))

#weights是权重，此文件中假设函数暂定为h(x)=1/1+exp(-w1*x1-w2*x2)
#未加入偏置b(就是w0*x0)，以便作图直观一些
#weights为2X10000的均匀序列，为要观察代价函数的变化

suppose = sigmoid(sampleDataMat.dot(weights))
print ('此处为假设函数输出的结果')
print (suppose)
print (np.shape(suppose))
print (suppose[0][0])

#随手检验一下假设函数的所有输出是否有问题
for i in range(len(suppose)):
    for j in range(len(suppose[i])):
        if suppose[i][j]>1 or suppose[i][j]<0:
            print ('有问题')

classLabelMat = np.tile(classLabelMat,(1,10000))
print ('将实际标签沿着1轴重复10000次，以便一百万个结果都能做差')
print (classLabelMat)

#下面注释掉的autoNorm函数为均值归一化
#加入了特征缩放会让图像变得更加平滑一点
#不加入均值归一化的话，图像震荡的会明显一些
#想看看加入效果的话，直接打开注释即可，无需进行其他操作
'''def autoNorm(dataSet):
    minVals = dataSet.min(0)
    maxVals = dataSet.max(0)
    ranges = maxVals - minVals
    normDataSet = np.zeros(np.shape(dataSet))
    m = dataSet.shape[0]
    normDataSet = dataSet - np.tile(minVals,(m,1))
    normDataSet = dataSet/np.tile(ranges,(m,1))
    return normDataSet,ranges,minVals

sampleDataMat,xxx,yyy = autoNorm(sampleDataMat)
'''

def costFunction(a,b):
    #下面为计算代价函数的具体代码过程
    #由于已知数据集有100个样本，遂以下进行重组数据结构的时候直接传入具体数字而非m了
    a = np.tile(a,(1,100))
    b = np.arange(100)
    b = np.tile(b,(100,1))
    b = b.transpose()
    b = b.reshape(1,10000)
    weights = np.vstack((a,b))
    costConsequence = sigmoid(sampleDataMat.dot(weights))-classLabelMat
    costConsequence = costConsequence**2
    costConsequence = costConsequence.sum(0)/200
    costConsequence = costConsequence.reshape((100,100))
    return costConsequence
 
#横纵坐标为0-99的均匀序列，X为权重w1，Y为权重w2    
X=np.arange(m)
Y=np.arange(m)

#绘制代价函数costFunction(w1,w2)的变化3D图
cost = costFunction(X,Y)
X,Y=np.meshgrid(X,Y)
fig=plt.figure()
ax = Axes3D(fig)
ax.plot_surface(X, Y, cost, rstride = 1, cstride = 1, cmap='rainbow')


    

最后是代码的执行结果，代价函数的图像（未进行特征缩放）：

从图中我们可以看出（尽管不是非常明显），在很多区域，确实存在多个波动，那么使用梯度下降等优化算法时，容易陷入局部最优解的确是个令人头疼的问题。