数据挖掘算法（三）--logistic回归

数据挖掘算法学习笔记汇总数据挖掘算法（一）–K近邻算法（KNN）数据挖掘算法（二）–决策树数据挖掘算法（三）–logistic回归在介绍logistic回归之前先复习几个基础知识点，有助于后面的理解。基本数学知识点1、对数似然函数若总体X为离散型，其概率分布列为 P(X=x)=p(x,θ) P(X=x)=p(x,\theta)其中θ\theta为未知参数。设 (X1,X2,...,

程序员学编程

1885人浏览 · 2017-09-26 19:18:36

程序员学编程 · 2017-09-26 19:18:36 发布

数据挖掘算法学习笔记汇总
数据挖掘算法（一）–K近邻算法（KNN）
数据挖掘算法（二）–决策树
 数据挖掘算法（三）–logistic回归

在介绍logistic回归之前先复习几个基础知识点，有助于后面的理解。

基本数学知识点

1、对数似然函数

若总体X为离散型，其概率分布列为

P (X = x) = p (x, θ)

θ <script type="math/tex" id="MathJax-Element-2">\theta</script>为未知参数。设

(X1,X2,...,Xn) <script type="math/tex" id="MathJax-Element-3"> (X_1,X_2,...,X_n) </script> 是取自总体样本容量为n的样本，则

(X1,X2,...,Xn) <script type="math/tex" id="MathJax-Element-4">(X_1,X_2,...,X_n)</script>的联合概率分布率为

\prod i = 1 n p (x i, θ)

(X1,X2,...,Xn) <script type="math/tex" id="MathJax-Element-6"> (X_1,X_2,...,X_n) </script>的一组观测值为

(x1,x2,...,xn) <script type="math/tex" id="MathJax-Element-7"> (x_1,x_2,...,x_n) </script>，易知样本

X1,X2,...,Xn <script type="math/tex" id="MathJax-Element-8"> X_1,X_2,...,X_n </script>取到观测值

x1,x2,...,xn <script type="math/tex" id="MathJax-Element-9"> x_1,x_2,...,x_n </script> 的概率为

L (θ) = L (x 1, x 2, . . ., x n; θ) = \prod i = 1 n p (x i, θ)

θ <script type="math/tex" id="MathJax-Element-11">\theta</script> 的取值而变化，它是

θ <script type="math/tex" id="MathJax-Element-12">\theta</script> 的函数，称

L(θ) <script type="math/tex" id="MathJax-Element-13">L(\theta)</script> 为样本的似然函数。但是由于来连乘的函数处理起来比较麻烦，所以对

L(θ) <script type="math/tex" id="MathJax-Element-14">L(\theta)</script> 取自然对数变成加法来处理要简单点。

l n L (θ) = \sum i = 1 n l n p (x i, θ)

2、logistic函数

logistic函数或logistic曲线是常见的“S”形（sigmoid curve ，S形曲线），方程式如下：

f (x) = L 1 + e - k ( x - x 0 )

e <script type="math/tex" id="MathJax-Element-17">e</script>自然对数
x0<script type="math/tex" id="MathJax-Element-18">x_0</script> S形中点的x值
L <script type="math/tex" id="MathJax-Element-19">L</script>曲线的最大值
k<script type="math/tex" id="MathJax-Element-20">k</script>曲线的陡度

上图是 L=1,k=1,x0=0 <script type="math/tex" id="MathJax-Element-21">L=1,k=1,x_0=0</script>时的图像
这里主要说明下这个函数的导数的性质，后面推导的时候会用到。
$f (x) = 1 1 + e - x = e x 1 + e x$ <script type="math/tex; mode=display" id="MathJax-Element-22">f(x)=\frac{1}{1+e^{-x}}=\frac{e^{x}}{1+e^{x}}</script>
$d d x f (x) = e x ( 1 + e x ) - e x e x ( 1 + e x ) 2$ <script type="math/tex; mode=display" id="MathJax-Element-23">\frac{d}{dx}f(x)=\frac{e^{x} (1+e^{x})-e^{x} e^{x}}{(1+e^{x})^2}</script>
$d d x f (x) = e x ( 1 + e x ) 2 = f (x) (1 - f (x))$ <script type="math/tex; mode=display" id="MathJax-Element-24">\frac{d}{dx}f(x)=\frac{e^{x}}{(1+e^{x})^2}=f(x)(1-f(x))</script>

logistic回归数学推导

先看一个简单的例子：
这里写图片描述
我们将平面上的点分为两类，中间的红色线条为边界。
预测类别 y=1 <script type="math/tex" id="MathJax-Element-25">y=1</script> 如果 −3+x1+x2≥0 <script type="math/tex" id="MathJax-Element-26">-3+x_1+x_2\geq0</script>预测类别 y=0 <script type="math/tex" id="MathJax-Element-27">y=0</script> 如果 −3+x1+x2<0 <script type="math/tex" id="MathJax-Element-28">-3+x_1+x_2 < 0</script>
此例子中

h θ (x) = g (θ 0 + θ 1 x 1 + θ 2 x 2)

对更多维的数据进行分类时，线性边界的情况，边界形式如下：

θ 1 x 1 + θ 2 x 2 + . . . + θ n x n = θ T x

h θ (x (i) ） = g (θ T x i) = 1 1 + e - θ T x i

hθ(x(i) <script type="math/tex" id="MathJax-Element-32">h_{\theta}(x^{(i)}</script>函数的值有特殊的含义，它表示结果取1的概率，因此对于输入x分类结果为类别1和类别0的概率分别为：

P (y = 1 | x; θ) = h θ (x (i)

P (y = 0 | x; θ) = 1 - h θ (x (i)

P (y | x ； θ) = (h θ (x)) y (1 - h θ (x)) 1 - y

L (θ) = \prod i = 1 m P (y (i) | x (i), θ)

l (θ) = l n L (θ)

= \sum i = 1 m (y (i) l n h θ (x (i)) + (1 - y (i)) l n (1 - h θ (x (i))))

l(θ) <script type="math/tex" id="MathJax-Element-39">l(\theta)</script>取最大值时的

θ <script type="math/tex" id="MathJax-Element-40">\theta</script>，利用梯度上升法求解，求得的

θ <script type="math/tex" id="MathJax-Element-41">\theta</script>就是要求的最佳参数。下面是利用梯度上升法求解过程。
求利用梯度上升法求解

l(θ) <script type="math/tex" id="MathJax-Element-42">l(\theta)</script>的最大值时，根据梯度上升法知道

θ <script type="math/tex" id="MathJax-Element-43">\theta</script>的更新公式如下：

θ j : = θ j + α \partial \partial θ j l (θ) (j = 0 . . . n)

l(θ) <script type="math/tex" id="MathJax-Element-45">l(\theta)</script>的偏导数：

\partial \partial θ j l (θ) = \sum i = 1 m ((y (i) 1 h θ ( x ( i ) ) \partial \partial θ j h θ (x (i)) - (1 - y (i)) 1 1 - h θ ( x ( i ) ) \partial \partial θ j h θ (x (i))

= \sum i = 1 m ((y (i) 1 g ( θ T x ( i ) ) - (1 - y (i)) 1 1 - g ( θ T x ( i ) )) \partial \partial θ j g (θ T x (i))

因为 g(θTxi) <script type="math/tex" id="MathJax-Element-48">g(\theta^Tx^{i})</script>是logistic函数

g (θ T x i) = 1 1 + e - θ T x i

l(θ) <script type="math/tex" id="MathJax-Element-50">l(\theta)</script>的偏导数转化

\partial \partial θ j l (θ) = \sum i = 1 m ((y (i) 1 g ( θ T x ( i ) ) - (1 - y (i)) 1 1 - g ( θ T x ( i ) )) g (θ T x (i)) (1 - g (θ T x (i))) \partial \partial θ j θ T x (i)

= \sum i = 1 m (y (i) (1 - g (θ T x (i))) - (1 - y (i)) g (θ T x (i))) x (i) j

= \sum i = 1 m (y (i) - g (θ T x (i))) x (i) j

= \sum i = 1 m (y (i) - h θ (x (i))) x (i) j

这样就得到了更新的过程

θ j : = θ j + α \sum i = 1 m (y (i) - h θ (x (i))) x (i) j (j = 0 . . . n)

python代码实现

本文代码运行环境：
python：3.5.1
pandas：0.19.2
其他环境可能有细微差别

# -*coding:utf-8*-
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import math

# 获取数据
data = pd.read_table("./logistic.txt", sep="\t", header=None)
dataMat = data.iloc[:, 0:-1]
labelMat = data.iloc[:, -1]


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

# 梯度上升算法
def gradAscent(dataMatrix, LabelsVector):
    n = dataMatrix.shape[1]
    alpha = 0.001
    maxCycles = 500
    thetas = np.ones((n, 1))
    for k in range(maxCycles):  # heavy on matrix operations
        h = sigmoid(dataMatrix * thetas)  # matrix mult
        error = LabelsVector.T - h  # vector subtraction
        thetas = thetas + alpha * dataMatrix.T * error  # matrix mult
    return thetas


def plotBestFit(thetas, data):
    """    
    :param thetas: type DataFrame , the thetas 
    :param data: type DtaFrame , all the data
    :return: 
    """
    X1 = data[data[3] == 0]
    X2 = data[data[3] == 1]
    fig = plt.figure()
    ax = fig.add_subplot(111)
    ax.scatter(X1[1], X1[2], s=30, c='red', marker='s')
    ax.scatter(X2[1], X2[2], s=30, c='green')
    x = np.arange(-3.0, 3.0, 0.1)
    y = (-thetas.iloc[0, 0] - thetas.iloc[1, 0] * x) / thetas.iloc[2, 0]
    ax.plot(x, y)
    plt.xlabel('X1')
    plt.ylabel('X2')
    plt.show()

thetas = gradAscent(np.mat(dataMat), np.mat(labelMat))
plotBestFit(pd.DataFrame(thetas), data)