机器学习实战（4）--logistic回归

1，sigmoid函数
对二分类的情形，类别C1的后验概率可以写成：

其中我们定义了：

且是logistic sigmoid函数，定义为
2，逻辑回归模型
考虑具有n个独立变量的向量x=(x1,x2,….,xn),设条件概率p(y=1|x)=p为根据观测相对于某件事发生的概率。由1的推算，逻辑回归模型可以表示为：

且。
g(x)一般表示为
其中w=(w0,w1,w2,….,wn),x=(1,x1,x2,….,xn).
3，最大似然函数
现在使用最大似然方法来确定逻辑回归模型中的参数w.
假设有n个观测样本，观测值分别y1,y2,….,yn.设pi=P(yi=1|xi)为给定条件下得到yi=1的概率,则p(yi=0|xi)=1-pi。故观测概率为P(yi)=pi^(yi)*(1-pi)^(1-yi).因为各个观测独立，所以他们的联合分布为各个边际分布的的乘积，即似然函数为：

（最大化似然函数的意义：我们的目的是找到一个参数序列w使得，目前观测到的数据y的概率最大。最大化概率函数就是认为已经观测到的数据具有最大的概率）。
直接求不好求，所以对似然函数取负对数，求对数似然函数的最小值。

4，最优化方法
在这里我们可以用两种方法来求参数值。
（1）梯度上升法
梯度的方向是函数值上升最快的方向，负梯度即为下降最快的方向。
Repeat Until convergence{

}
其中，w是参数，t表示类别1或0，y即为根据参数以及x计算出来的值。梯度上升法每一次都用到了全部的样本数据。
（2）随机梯度法
Loop{
for j to n：
{

}

}
随机梯度上升法每次只用了样本中的一个值。
对于梯度下降以及随机梯度下降只要将a前面的改为负号即可。
5，分类
通过计算出参数值w，这样每次新来的一个值x，只要计算
和，判断哪个值大于0.5，即x属于该类。

6，Python实现

#-*- coding:gb2132 -*-from numpy import *def loadDataSet():    dataMat = []; labelMat = []    fr = open('testSet.txt')    for line in fr.readlines():        lineArr = line.strip().split()        #x=(1,x1,x2,...,xn)        dataMat.append([1.0, float(lineArr[0]), float(lineArr[1])])        labelMat.append(int(lineArr[2]))    return dataMat,labelMatdef sigmoid(inX):    return 1.0/(1+exp(-inX))def gradAscent(dataMatIn, classLabels):    dataMatrix = mat(dataMatIn)             #convert to NumPy matrix    labelMat = mat(classLabels).transpose() #convert to NumPy matrix    m,n = shape(dataMatrix)    alpha = 0.001    maxCycles = 500    #初始化权值为1    weights = ones((n,1))    for k in range(maxCycles):              #heavy on matrix operations        h = sigmoid(dataMatrix*weights)     #matrix mult        error = (labelMat - h)              #vector subtraction        #梯度下降法        weights = weights + alpha * dataMatrix.transpose()* error #matrix mult    return weights#画决策线def plotBestFit(weights):    import matplotlib.pyplot as plt    dataMat,labelMat=loadDataSet()    dataArr = array(dataMat)    n = shape(dataArr)[0]     xcord1 = []; ycord1 = []    xcord2 = []; ycord2 = []    for i in range(n):        if int(labelMat[i])== 1:            xcord1.append(dataArr[i,1]); ycord1.append(dataArr[i,2])        else:            xcord2.append(dataArr[i,1]); ycord2.append(dataArr[i,2])    fig = plt.figure()    ax = fig.add_subplot(111)    ax.scatter(xcord1, ycord1, s=30, c='red', marker='s')    ax.scatter(xcord2, ycord2, s=30, c='green')    x = arange(-3.0, 3.0, 0.1)    y = (-weights[0]-weights[1]*x)/weights[2]    ax.plot(x, y)    plt.xlabel('X1'); plt.ylabel('X2');    plt.show()#随机梯度法def stocGradAscent0(dataMatrix, classLabels):    m,n = shape(dataMatrix)    alpha = 0.01    weights = ones(n)   #initialize to all ones    for i in range(m):        h = sigmoid(sum(dataMatrix[i]*weights))        error = classLabels[i] - h        weights = weights + alpha * error * dataMatrix[i]    return weights#随机梯度法，不过a值在学习过程中不断改变def stocGradAscent1(dataMatrix, classLabels, numIter=150):    m,n = shape(dataMatrix)    weights = ones(n)   #initialize to all ones    for j in range(numIter):        dataIndex = range(m)        for i in range(m):            alpha = 4/(1.0+j+i)+0.0001    #apha decreases with iteration, does not             #随机取值进行训练            randIndex = int(random.uniform(0,len(dataIndex)))#go to 0 because of the constant            h = sigmoid(sum(dataMatrix[randIndex]*weights))            error = classLabels[randIndex] - h            weights = weights + alpha * error * dataMatrix[randIndex]            del(dataIndex[randIndex])    return weightsdef classifyVector(inX, weights):    prob = sigmoid(sum(inX*weights))    if prob > 0.5: return 1.0    else: return 0.0def colicTest():    frTrain = open('horseColicTraining.txt'); frTest = open('horseColicTest.txt')    trainingSet = []; trainingLabels = []    for line in frTrain.readlines():        currLine = line.strip().split('\t')        lineArr =[]        for i in range(21):            lineArr.append(float(currLine[i]))        trainingSet.append(lineArr)        trainingLabels.append(float(currLine[21]))    trainWeights = stocGradAscent1(array(trainingSet), trainingLabels, 1000)    errorCount = 0; numTestVec = 0.0    for line in frTest.readlines():        numTestVec += 1.0        currLine = line.strip().split('\t')        lineArr =[]        for i in range(21):            lineArr.append(float(currLine[i]))        if int(classifyVector(array(lineArr), trainWeights))!= int(currLine[21]):            errorCount += 1    errorRate = (float(errorCount)/numTestVec)    print "the error rate of this test is: %f" % errorRate    return errorRate#运行10次，测试平均误差def multiTest():    numTests = 10; errorSum=0.0    for k in range(numTests):        errorSum += colicTest()    print "after %d iterations the average error rate is: %f" % (numTests, errorSum/float(numTests))

horseColicTraining.txt：http://pan.baidu.com/s/1bogV7aR
horseColicTest.txt:http://pan.baidu.com/s/1NuF0a
testSte.txt:http://pan.baidu.com/s/1eQRWjOA