碎碎念

好久不更新，四月重新开始！

台大李宏毅机器学习作业——分类算法

Posted on Edited on In
Symbols count in article: 2.7k Reading time ≈ 2 mins.

好久没更，李宏毅第二次课程实操

第二次作业主要做的是分类算法。分类算法分为Generative和Discriminate两种。Generative常见的比如朴素贝叶斯，Discriminator就多了，逻辑回归就是。这里实现两个不同的分类算法。

首先是probabilistic generative model。概率生成模型其实就是根据数据分布的情况，拟合一种分布，计算条件概率，判断样本的类别。基本公式就是\(P(C_i | x) = \frac{P(x | C_i) P(C_i)}{\sum(P(x | C_i) P(C_i))}\)。那需要用到概率分布的地方就是\(P(x | C_i)\)这里，我们需要假设一个分布来拟合这个概率。一般而言，我们用正态分布比较多。之所以用正态分布多，可以回顾一下大数定律和中心极限定理。当数据量到了一定程度，就会呈现出正态分布。

大概的原理可以看这篇博客，不是一个很复杂的模型。这里用到的数据就是作业2的示例数据，是一个二分类的数据。因为模型的公式已经很明显了，这里就直接将公式转化为代码：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
def gen_model(X_test, X_train, threshold=0.5):
    """Generative model"""
    D = X_train.shape[1] / 2
    a = X_train.loc[X_train['y'] == 0, :].drop('y', axis=1)
    mu1 = a.mean()
    tmp1 = pd.DataFrame(X_test - mu1)
    b = X_train.loc[X_train['y'] != 0, :].drop('y', axis=1)
    mu2 = b.mean()
    tmp2 = pd.DataFrame(X_test - mu2)
    # 这里的协方差矩阵是两个类别的协方差矩阵的加权平均。通常来说用加权平均的效果比较好。
    sigma = a.shape[0] / len(X) * np.cov(a.T, bias=True) + b.shape[0] / len(X_train) * np.cov(b.T, bias=True)
    sigma_det = 1 / np.sqrt(np.linalg.det(sigma))
    # 这里用伪逆，实际上不用伪逆问题也不大
    sigma_pinv = np.linalg.pinv(sigma)
    # 这里要注意，我直接算的是矩阵，而不是一行一行的算。那么直接计算矩阵的时候其实就只要拿对角线元素就可以了。为什么直接取对角线元素，可以回顾一下MIT的线性代数课程。
    array1 = np.diag(np.dot(tmp1, sigma_pinv).dot(tmp1.T))
    array2 = np.diag(np.dot(tmp2, sigma_pinv).dot(tmp2.T))
    gauss1 = np.power(2 * np.pi, -D) * sigma_det * np.exp(-0.5 * array1)
    gauss2 = np.power(2 * np.pi, -D) * sigma_det * np.exp(-0.5 * array2)
    prob1 = gauss1 * a.shape[0] / len(X)
    prob2 = gauss2 * b.shape[0] / len(X)
    prob = prob1 / (prob1 + prob2)
    prob[np.isnan(prob)] = 0
    # 这里用了一个阈值，其实也可以按照每个类别的先验概率来分类。个人以为比较合理，但是实际上未必效果最好。五五开也是一种常见的方法。
    prob[prob >= threshold] = 0
    prob[prob != 0] = 1
    return prob

那实际上这样概率生成模型就结束了。从某种程度而言，概率生成模型是不存在训练测试的。

这里顺带提一下用的几个包：

1
2
3
4
5
6
import numpy as np
import pandas as pd
from sklearn.metrics import classification_report
from sklearn.metrics import confusion_matrix
from sklearn.metrics import accuracy_score
import time

其实sklearn的几个metric可以自己写，但是这里偷懒了。原则上还是自己实现核心算法。

相对而言，逻辑回归就比较简单了，大概过程跟线性回归是一样的，用梯度下降去找\(w\)。原理看这篇博客。所以实现也很简单：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
def sigma(x):
    return 1 / (1 + np.exp(-x))


def minibatch_sgd(X, y_true, w, eta=0.1, epoch=10, batch_size=100):
    rounds = 0
    while rounds < epoch:
        start = time.time()
        start_index = 0
        sum_error = 0
        while start_index < (X.shape[0] - 1 - batch_size):
            end_index = start_index + batch_size
            if end_index > X.shape[0]:
                end_index = X.shape[0] - 1
            X_batch = X.iloc[start_index:end_index, :]
            error = sigma(np.dot(X_batch, w)) - y_true[start_index:end_index].values
            sum_error += sum(error)
            w -= eta * np.dot(X_batch.T, error)
            start_index += batch_size
        rounds += 1
        end = time.time()
        accuracy = accuracy_score(sigma(np.dot(X, w)), y_true)
        print('epoch: ' + str(rounds) + '  time: '+ str(round(end - start, 2)) + 's' + '   error: ' + str(sum_error) + '  accuracy: ' + str(accuracy))
    return(w)

这里用的是minibatch_sgd，速度会快一点。但是不知道为什么我这边跑的每一个epoch的loss都会乱荡，有点尴尬，后面再仔细debug。