Logistic回归

一、logistic回归

在吴恩达教授的深度学习中的logistic回归是通过一个概率模型来对样本进行分类。比如，输入特征$x$和相应的标签$y$（$y$是$0$或$1$，表示两个类别）。我们希望找到一个参数化的模型，可以最好地描述$y$在给定$x$的条件下的分布。
在logistic中，用sigmoid函数来建模概率分布。

$$\sigma(z) = \frac{1}{1 + e^{-z}}$$

二、logistic回归实现分类猫咪图片

$X$、$Y$

我是这样理解的，就拿吴恩达教授课上的判断是不是猫的图片来举例吧。如果要实现对图片分类，首先输入是肯定是图片，然后输出就是这张图片 $是/不是$ 猫咪。
对于输入，一张彩色的猫咪图片的颜色信息可以分为$RGB$三个通道，如果是对于一张$64*64$像素的彩色猫咪图片，我们可以把它的信息整理为一个$(64,64,3)$的矩阵$x$。如果是对于很多张图片，我们可以把它整成一个$(m,64,64,3)$的矩形，我们就把这个矩阵叫$X$吧。

$$X = [x^{(1)} , x^{(2)} , x^{(3)} , ...... ,x^{(i)}]$$

对于输出就很简单了，判断一个图片是不是猫咪就行了，我们用0来表示不是猫咪，1来表示是猫咪。所以可以很容易就可以知道$y = 0 | y = 1$。对于一系列图片的输入，我们就用一个$(m,1)$的矩阵$Y$来表示每张图片对应的输出就行了

$$Y = [y^{(1)} , y^{(2)} , y^{(3)} , ...... ,y^{(i)}]$$ $$对于任意y^{(j)}，y^{(j)} = 1 或y^{(j)} = 0$$

然后接下来就是用训练集训练模型了，其实所谓的训练模型就是拿一组已经确定好是不是猫咪的图片（Y是确定的）拿去找X和Y之间的关系，然后得到相应的参数。举个帮助理解的例子，对于一个式子$y = ax$，我们知道了$y$和$x$就可以求对应的$a$。
然后是测试模型，测试就相当于实操了，给你一组图片让你模型去判断是不是真的猫咪。这时候相当于只知道X和训练出来的那些参数，去预测这个是不是Y。相当于知道了$x$和$a$只要求$y$就行了。
不过logistic回归并没有$y = ax$这么简单直接。他是其实知道了$X$和训练出来的参数之后并不是直接求出$Y$，而是给出一个$Y$可能是1的概率。然后我们可以给他加一个阈值，这个概率如果大于多少多少，我们就可以直接判断这张图片就是一只猫。
$Y$(这个Y代表的是随机变量，不是矩阵$Y$)其实是一个0-1分布：

$Y$	$0$	$1$
$P$	$1-\hat{y}$	$\hat{y}$
$\hat{y}$代表的就是$y$是$1$的概率

$$P\{y^{(i)}|x^{(i)}\} = \hat{y}^y(1-\hat{y})^y,y = 0,1,0<\hat{y}<1.$$

线性组合和sigmoid函数

在logistic回归中，我们用一个线性公式$z = w^TX + b$来对输入样本X进行线性组合，然后将得到的实数z映射到0到1之间的概率值(当然，对于单个特征z是实数，但是有m个特征的z就是一个(1，m)的矩阵)。我们就用这个sigmoid函数来映射：

$$P(y = 1 | x; w, b) = \hat{y} = A = \sigma(z) = \frac{1}{1 + e^{-z}}$$

其中

$$z = w^TX + b$$

w代表权重,是一个维度为(m,1)的向量。b代表偏移量，是一个实数。样本X是一个维度为$(64×64×3,m)$的矩阵。

一开始X是一个(m,64,64,3)的矩阵，不过这不方便我们线性组合，所以我们就把X降维成(64×64×3,m)，其中降维方法就是先行后列拼接起来而已。
最后我们训练其实就是为了训练出最合适的w和b，然后根据这个w和b去预测Y。

最大似然函数和损失函数

所以我们就可以用最大似然函数来通过观测到的数据来估计参数值，使得观测数据出现的概率最大。在这里其实就是找到一组$w$和$b$，然后在给定的训练样本$(X,Y)$下，观测到的$y^{(i)}$出现的概率最大

$$L(w,b)=\displaystyle \prod^{m}_{i = 1}{P(y(i) | x(i); w, b)}$$

由于$\ln{x}$是$x$的单调递增函数，所以$\ln{L(w,b)}$和$L(w,b)$有相同的最大值点，所以我们直接用$\ln{L(w,b)}$，计算起来会比较方便。同时我们可以适当的放缩一下。

$$\cal L(w,b) = \it \frac{1}{m} \ln{L(w,b)} = \displaystyle \frac{1}{m}\sum^{m}_{i = 1}{\ln{P(y(i) | x(i);w,b)}}$$ $$=\frac{1}{m} \displaystyle \sum^{m}_{i = 1}{y^{(i)}\ln{\hat{y}^{(i)}}} + (1 - y^{(i)}) \ln{(1 - \hat{y}^{(i)})}$$

不过这个图像是凸的，为了符合直观，我们就把这个式子加上一个负号，这样图像就变成凹的了,求最大值就变成了求最小值。
我们把加负号的式子叫做损失函数，用$J(w,b)$来表示

$$J(w,b) = -\cal L(w,b)$$

梯度下降求取最优参数

在学习最大释然函数的时候，都是直接令偏导为0，然后求出最优的参数，但是在Logistic回归中，我们是没有办法求出最优的参数的，所以我们就利用梯度下降一步一步的来接近最优的参数。

$$w = w - \alpha \frac{\partial(J(w,b))}{\partial w}$$ $$b = b - \alpha \frac{\partial(J(w,b))}{\partial b}$$

其中

$$\frac{\partial(J(w,b))}{\partial w} = \frac{\partial(J(w,b))}{\partial \hat{y}}\frac{\partial\hat{y}}{\partial z} \frac{\partial z}{\partial w} = \frac{1}{m}X(\hat{y} - y)$$ $$\frac{\partial(J(w,b))}{\partial b} = \frac{\partial(J(w,b))}{\partial \hat{y}}\frac{\partial\hat{y}}{\partial z} \frac{\partial z}{\partial b} = \frac{1}{m}(\hat{y} - y)$$

因为J(w,b)的图像是凹的，我们的w要达到最小值就要一步一步的接近最小值的位置，我们减去dw可以让w根据图像的斜率动态的变化，最后达到最小值的位置(想象一下，如果w在最小值左边，这时候dw<0,那么w = w - dw就会让w变大，使得w更加靠近最小值。w在右边同理)。这里的$\alpha$是学习率，代表每一次下降走的距离，如果$\alpha$太大，那么步子迈太大可能就直接跨过最小值了。而$\alpha$如果太小，那么可能就需要迈特别多步，才能到达最小值，就很浪费时间。所以$\alpha$的选择比较重要。

预测

最后，训练出来的$w$和$b$就可以用来预测了，对于一个新的样本$X$，我们就可以根据$\hat{y} = \sigma(w^T X + b)$来获取$y = 1$的概率，然后根据我们自己设定的阈值来判断这个样本是不是$1$

源代码

import numpy as np
import matplotlib.pyplot as plt
import h5py
from lr_utils import load_dataset


# 导入数据集

train_set_x_orig , train_set_y , test_set_x_orig , test_set_y , classes = load_dataset()

# 测试一下,看看训练集里的图
index = 50
plt.imshow(train_set_x_orig[index])
plt.show()

print("y=" + str(train_set_y[:,index]) + ", it's a " + classes[np.squeeze(train_set_y[:,index])].decode("utf-8") + "' picture")


m_train = train_set_y.shape[1] #训练集里图片的数量。
m_test = test_set_y.shape[1] #测试集里图片的数量。
num_px = train_set_x_orig.shape[1] #训练、测试集里面的图片的宽度和高度（均为64x64）。

# 我们要把xxx,64,64,3构造为64*64*3,xxxx
train_set_x_flatten  = train_set_x_orig.reshape(train_set_x_orig.shape[0],-1).T  #把四维降到两维，行数是209（图片数量209），列数-1让系统自己判断，然后转置变成 64*64*3行，209列
test_set_x_flatten = test_set_x_orig.reshape(test_set_x_orig.shape[0],-1).T   # 同理
print(train_set_x_flatten.shape)
print(test_set_x_flatten.shape)

# 然后对数据集进行标准化，因为RGB最大值也就255，所以把它们直接除以255可以保证数据都在0到1之间
train_set_x = train_set_x_flatten / 255
test_set_x =test_set_x_flatten / 255

# 构建sigmoid()
def sigmoid(z):
    "sigmoid(z) = 1 /(1 + e^(-z))"
    s = 1 / (1 + np.exp(-z))
    return s

# 然后初始化参数w和b
def init_with_zeros(dim):
    """
        初始化w为维度为(dim,1)0向量,初始化b为0
    """
    w = np.zeros(shape=(dim,1))
    b = 0

    # 同时确保我们的数据是正确的
    assert(w.shape == (dim,1))
    assert(isinstance(b, float) or isinstance(b, int)) # b得是int或者float才行

    return (w, b)


# 计算成本函数
def propagate(w, b, X, Y):
    """
        w 代表权重
        b 代表偏差
        X 数据集 维度(64*64*3, 训练数量)
        Y 输出，如果是猫就1，不是猫就0 维度(1, 训练数量)
    """

    m = X.shape[1] #获取数据集个数m

    # 计算y_hat
    z = np.dot(w.T,X) + b
    A = sigmoid(z)
    #成本函数 J
    cost = (-1 / m) * np.sum(Y * np.log(A) + (1 - Y) * np.log(1 - A))

    # J对w和b求偏导
    dw = (1 / m) * np.dot(X, (A - Y).T)
    db = (1 / m) * np.sum(A - Y)

    # 确保数据是正确的
    assert(dw.shape == w.shape)
    assert(db.dtype == float)
    cost = np.squeeze(cost)
    assert(cost.shape == ())

    # 用字典把dw和db保存起来
    grads = {
            "dw": dw,
            "db": db
            }

    return (grads, cost)

# 我们目标是通过最小化成本函数J来学习w和b。比如对于参数w，更新规则是w = w - α*dw，其中α是学习率。可以发现，根据曲线的递增递减对应的导数的正负，可以动态的调整w，让w越来越接近目标值
# 通过已知的X和Y来学习w和b，我们需要做的是最小化成本函数J，对于参数w，更新规则是w = w - α*dw，其中α是学习率。可以发现，根据曲线的递增递减对应的导数的正负，可以动态的调整w，让w越来越接近目标值
def optimize(w, b, X, Y, num_iterations, learning_rate, print_cost = False):
    """
        w,b,X,Y含义同上
        num_iterations 优化循环的迭代次数
        learning_rate  梯度下降更新规则的学习率
        print_cost  打印损失值
    """

    costs = []

    for i in range(num_iterations):
        # 梯度下降 最精华的部分
        grads, cost = propagate(w, b, X, Y)

        dw = grads["dw"]
        db = grads["db"]

        w = w - learning_rate * dw
        b = b - learning_rate * db

        # 记录成本
        if i % 100 == 0:
            costs.append(cost)  #c++里的push，每迭代100次记录一下，记录的数据放在costs里面
        # 打印成本数据
        if (print_cost) and (i % 100 == 0):
            print("迭代的次数: %i ， 误差值： %f" % (i, cost))


        params = {
                "w" : w,
                "b" : b
                }

        grads = {
                "dw": dw,
                "db": db
                }
    return (params, grads, costs)


# 然后用学习的w和b来预测数据集X的标签是0还是1,相当于知道了w和b和X，就可以去预测Y
def predict(w, b, X):
    """
        然后就是用逻辑回归参数来预测标签是0还是1了
    """

    m = X.shape[1] # 图片数量
    Y_prediction = np.zeros((1,m)) # 初始化
    print("reshape前：" + str(w.shape))
    w = w.reshape(X.shape[0], 1) # 让w的行数和X的参数数量一样，w(i)和x(i)一一对应
    print("reshape后：" + str(w.shape))
    # 预计猫在图片中出现的概率
    z = np.dot(w.T, X) + b
    A = sigmoid(z)
    for i in range(A.shape[1]):
        # 如果概率大于0.5就认为是猫
        if A[0, i] > 0.5:
            Y_prediction[0, i] = 1
        else:
            Y_prediction[0, i] = 0

    # 检查预测的Y
    assert(Y_prediction.shape == (1,m))

    return Y_prediction



# 最后是吧上面的函数都整合一下
def model(X_train, Y_train, X_test, Y_test, num_iterations=2000, learning_rate=0.5, print_cost=False):
    """
    通过调用之前实现的函数来构建逻辑回归模型

    参数：
        X_train  - numpy的数组,维度为（num_px * num_px * 3，m_train）的训练集
        Y_train  - numpy的数组,维度为（1，m_train）（矢量）的训练标签集
        X_test   - numpy的数组,维度为（num_px * num_px * 3，m_test）的测试集
        Y_test   - numpy的数组,维度为（1，m_test）的（向量）的测试标签集
        num_iterations  - 表示用于优化参数的迭代次数的超参数
        learning_rate  - 表示optimize（）更新规则中使用的学习速率的超参数
        print_cost  - 设置为true以每100次迭代打印成本

    返回：
        d  - 包含有关模型信息的字典。
    """

    # 初始化w和b
    w, b = init_with_zeros(X_train.shape[0])

    # 其实最关键的就是grads里的w和b，导入数据集，告诉他学习率和迭代次数，就可以得到w和b
    parameters, grads, costs = optimize(w, b, X_train, Y_train, num_iterations, learning_rate, print_cost)

    # 从字典“参数”中检索更新后的参数w和b
    w, b = parameters["w"], parameters["b"]

    # 预测测试/训练集的例子
    # 然后通过学习到的w和b来预测X_train和X_test
    Y_prediction_test = predict(w, b, X_test)
    Y_prediction_train = predict(w, b, X_train)

    # 打印训练后的准确性
    print("训练集准确性：", format(100 - np.mean(np.abs(Y_prediction_train - Y_train)) * 100), "%")
    print("测试集准确性：", format(100 - np.mean(np.abs(Y_prediction_test - Y_test)) * 100), "%")

    d = {
        "costs": costs,
        "Y_prediction_test": Y_prediction_test,
        "Y_prediciton_train": Y_prediction_train,
        "w": w,
        "b": b,
        "learning_rate": learning_rate,
        "num_iterations": num_iterations}
    return d


d = model(train_set_x, train_set_y, test_set_x, test_set_y, num_iterations = 2000, learning_rate = 0.005, print_cost = True)



#绘制图
costs = np.squeeze(d['costs'])  # squeeze会消去一个轴，比如对于数组[1]，用了squeeze就会变成纯数值的1了，对于数值来说，他的shape是空的
plt.plot(costs)
plt.ylabel('cost')
plt.xlabel('iterations (per hundreds)')
plt.title("Learning rate =" + str(d["learning_rate"]))
plt.show()

哎，这样学累死了，会用就行了。

[1]

【中文】【吴恩达课后编程作业】Course 1 - 神经网络和深度学习 - 第二周作业