可以学习课程内容李宏毅机器学习特训营

项目描述

二元分类是机器学习中最基础的问题之一，在这份教学中，你将学会如何实作一个线性二元分类器，来根据人们的个人资料，判断其年收入是否高于 50,000 美元。

我们将以两种方法: logistic regression 与 generative model，来达成以上目的，你可以尝试了解、分析两者的设计理念及差别。 实现二分类任务：

个人收入是否超过50000元？

数据集介绍

这个资料集是由UCI Machine Learning Repository 的Census-Income (KDD) Data Set 经过一些处理而得来。为了方便训练，我们移除了一些不必要的资讯，并且稍微平衡了正负两种标记的比例。事实上在训练过程中，只有 X_train、Y_train 和 X_test 这三个经过处理的档案会被使用到，train.csv 和 test.csv 这两个原始资料档则可以提供你一些额外的资讯。

已经去除不必要的属性。
已经平衡正标和负标数据之间的比例。

特征格式

train.csv，test_no_label.csv。
基于文本的原始数据
去掉不必要的属性，平衡正负比例。
X_train, Y_train, X_test(测试)
train.csv中的离散特征=>在X_train中onehot编码(学历、状态…)
train.csv中的连续特征 => 在X_train中保持不变(年龄、资本损失…)。
X_train, X_test : 每一行包含一个510-dim的特征，代表一个样本。
Y_train: label = 0 表示 “<=50K” 、 label = 1 表示 \" >50K \" 。

项目要求

请动手编写 gradient descent 实现 logistic regression
请动手实现概率生成模型。

单个代码块运行时长应低于五分钟。
禁止使用任何开源的代码(例如，你在GitHub上找到的决策树的实现)。

logistic regression

先导入数据

import numpy as np

X_train_path = \'work/data/X_train\'
Y_train_path = \'work/data/Y_train\'
X_test_path = \'work/data/X_test\'
with open(X_train_path) as f:
    next(f)
    X_train = np.array([line.strip(\'\\n\').split(\',\')[1:] for line in f],dtype=float)
with open(Y_train_path) as f:
    next(f)
    Y_train = np.array([line.strip(\'\\n\').split(\',\')[1] for line in f], dtype=float)
with open(X_test_path) as f:
    next(f)
    X_test = np.array([line.strip(\'\\n\').split(\',\')[1:] for line in f], dtype=float)

print(X_train.shape)
print(Y_train.shape)
print(Y_train[0:3]) 

先把X_train 和 Y_train 组合成一个矩阵，构成一个数据集data，最后一列是真实标签，这一个操作用于方便后面的随机打乱操作

Y_train=Y_train.reshape(-1,1)
data = np.concatenate((X_train,Y_train),axis=1)
print(data.shape)

写一个随机打乱数据的函数，下面函数的就是按行随机打乱

#按行随机打乱训练数据集
def shuffle(x):
    return np.random.permutation(x)

对数据进行标准化 z-score normalization 的想法

在训练之前，我们都会输入数据进行标准化处理，公式就是 （x - μ）/σ，我这里不谈为什么要进行标准化。 网上的有些作业方法里面，很多是首先对数据直接进行标准化，就比如上面有了X_train，就直接运行了下面的代码，对所有数据标准化之后，在这个新获得的X_train的基础上进行训练集和验证集的划分，其实这里的划分我觉得是有一点不合理的。

X_mean = np.mean(X_train ,axis=0)
X_mean = X_mean.reshape(1,510)
X_std  = np.std(X_train, axis=0).reshape(1, -1)

X_train = (X_train-X_mean)/(X_std+0.00000001) #避免方差为0，当方差为零时，说明这一列的数值是一样的，那就把这一列都变为0
X_train

划分出验证集的目的就是评价一下当前模型的好坏，以便于根据情况调整超参数，让模型训练的更好，那么验证集就应该尽可能的与训练集没有依赖关系。那么合理的做法应该是先对原始的X_train 进行随机划分——分为训练集和验证集，然后分别对训练集和验证集进行标准化。

这里说一下我的原因：对于训练出来的模型，以后的输入数据是完全未知的，最后得到的准确率可能没有想象中的高。如果一次性标准化整个数据集，就忽略了新数据输入会超出标准化范围的可能性。

如果先对X_train标准化，然后划分训练集和验证集，那此时的验证集对于模型来说是已经“见到过了”，或者说这个验证集对模型来说并不是未知的，因为验证集的一些数据信息已经平均在训练集中了（说法可能不准确），如果这样的验证集去评价模型，可能在训练过程中，模型的表型还可以，但当模型去面对一个未知的数据时，可能表现就差了。

合理的做法是要确保验证数据是“未知的”，不要通过任何方式“偷看试卷答案”！先对原始的X_train
进行随机划分——分为训练集和验证集，然后分别对训练集和验证集进行标准化！

#对数据集data进行随机打乱
data_shuffle = shuffle(data)
#进行训练集和验证集划分
train_len = int(len(data)*0.9)  #这里训练集占90%
train_set = data_shuffle[0:train_len,:]
vali_set = data_shuffle[train_len:,:]
print(train_set.shape)
print(vali_set.shape)

#下面对训练集和验证集分别进行标准化
train_set_mean = np.mean(train_set[:,0:-1] ,axis=0)
train_set_mean = train_set_mean.reshape(1,510)
print(train_set_mean.shape)
train_set_std  = np.std(train_set[:,0:-1], axis=0).reshape(1, -1)
train_set[:,0:-1]  = (train_set[:,0:-1]-train_set_mean)/(train_set_std+0.00000001) #避免方差为0。当方差为零时，说明这一列的数值是一样的，那就把这一列都变为0

vali_set[:,0:-1] = (vali_set[:,0:-1] - train_set_mean)/(train_set_std+0.00000001)  #对验证集进行标准化时，使用训练集的均值和方差

print(train_set.shape)
print(vali_set.shape)
print(train_set[0,:])

做logistic regression，需要Sigmoid Function，可以写一个Sigmoid Function函数。

然后需要一个损失函数，cross entropy

明确了需要的函数后，就可以开始写函数了。

#Sigmoid Function函数
def sigmoidFun(z):
    ##########################################################
    ##  写Sigmoid 函数可能会出现overflow encountered 的警告
    ##  例如：import numpy as np                            
    ##        x = np.array([-10000,2,3,5,6,7,8,900])
    ##        1.0/(1+np.exp(-x))                            
    ##  运行就会有警告/opt/conda/envs/python35-paddle120-env/lib/python3.7/site-packages/ipykernel_launcher.py:3: RuntimeWarning: overflow encountered in exp
    ##   This is separate from the ipykernel package so we can avoid doing imports until
    ##  这是因为numpy数组x中可能有绝对值比较大的负数，这样传给sigmoid函数时，分母np.exp(-x)会非常大，导致np.exp(-x)溢出
    ##  在运算之前，可以根据情况先对数组限制一下最小最大值。
    return 1/(1.0+np.exp(-np.clip(z,-10,10)))

#损失函数cross entropy
def crossEntropy(y_pre,y_true):
    ## y_pre 为预测值
    ## y_true 为真实值    
    loss = -np.sum(y_true * np.log(y_pre) + (1 - y_true) * np.log(1 - y_pre))
    return loss

#计算预测准确率
def accuracy(Y_pred, Y_label):
    acc = 1 - np.mean(np.abs(Y_pred - Y_label))
    return acc

数据量很大时，要考虑用 Mini_Batch Gradient Descent（最小批量梯度下降）

假如加入全部训练数据有1000条，因为同时放进去训练计算量太大，那就每次只选择10条数据用于训练

下面的训练代码里，使用了Adaptive Learning Rates

解释一下for t in range(iter_time):中的代码

iter_time是训练的总轮数，每次在新一轮训练开始之前，都随机打乱一下训练集，因为训练集比较大，就用Mini_Batch Gradient Descent， batch就是一次送去训练的训练集大小。

batch = 1000
batch_times = int(train_len/batch)
lr = 0.1 #学习率
iter_time = 500  #训练的总轮数
epochs = 2
times = 1 #训练次数，学习率随着训练次数的增加而减小 
dim = 510 + 1 # 510是train_set_x的列数，1是bias（b）的维数
w = np.zeros([dim, 1]) #初始化 这里面的w已经包括了b   也可以分开初始化，w=np.zeros([510, 1]) b = np.zeros([1, 1])

for t in range(iter_time):

    train_set = shuffle(train_set)  # 随机打乱训练数据
    train_set_x = train_set[:,0:-1] #获取特征数据
    train_set_x = np.concatenate((train_set_x,np.ones([train_len,1])),axis=1)
    train_set_y = train_set[:,-1]   #获取标签值

    vali_set = shuffle(vali_set)    #随机打乱验证数据
    vali_set_x = vali_set[:,0:-1]  #获取验证集的特征数据
    vali_set_x = np.concatenate((vali_set_x,np.ones([len(data)-train_len,1])),axis=1) 
    vali_set_y = vali_set[:,-1]#获取验证集的标签数据
    vali_set_y = vali_set_y.reshape(-1,1)

    y = np.zeros([batch, 1]) #用于存储训练时的预测值
    y_hat = np.zeros([batch, 1]) #用于存储小批次的标签值

    for b in range(batch_times):
        x= train_set_x[batch*b:batch*(b+1),:]
        x = x.reshape(batch,-1)
        y_hat = train_set_y[b*batch:(b+1)*batch]
        y_hat = y_hat.reshape(batch,1)
        y = sigmoidFun(np.dot(x,w))  #用训练集 train_set_x 预测的y
        err = y - y_hat
        gradient = np.dot(x.transpose(),err) #求梯度
        w = w - lr/np.sqrt(times) * gradient  # 更新参数，Adaptive Learning Rates
        times = times+1

    if(t%100==0):  #每100轮打印一次
        y_predict = sigmoidFun(np.dot(vali_set_x,w))
        loss_vali  = crossEntropy(y_predict,vali_set_y)/(len(data)-train_len) #计算验证集交叉熵
        acc_vali = accuracy(np.round(y_predict),vali_set_y) #计算验证集准确率
        loos = crossEntropy(y,y_hat)/train_len #计算训练集交叉熵
        acc = accuracy(np.round(y),y_hat) #计算训练集准确率
        

        print(str(t)+\"/\" +str(iter_time) + \" 训练集交叉熵：\"+str(loos))
        print(str(t)+\"/\" +str(iter_time) + \" 训练集准确率：\"+str(acc))
        print(str(t)+\"/\" +str(iter_time) + \" 验证集交叉熵：\"+str(loss_vali))
        print(str(t)+\"/\" +str(iter_time) + \" 验证集准确率：\"+str(acc_vali))

训练结果：

测试

#对测试集进行标准化

test_set  = (X_test-train_set_mean)/(train_set_std+0.00000001) #还是用训练集的均值和方差

#这里对测试集增加一列
test_set = np.concatenate((test_set,np.ones([len(test_set),1])),axis=1)

predict = np.round(sigmoidFun(np.dot(test_set,w)))

Porbabilistic generative model

导入数据

import numpy as np

X_train_path = \'work/data/X_train\'
Y_train_path = \'work/data/Y_train\'
X_test_path = \'work/data/X_test\'
with open(X_train_path) as f:
    next(f)
    X_train = np.array([line.strip(\'\\n\').split(\',\')[1:] for line in f], dtype = float)
with open(Y_train_path) as f:
    next(f)
    Y_train = np.array([line.strip(\'\\n\').split(\',\')[1] for line in f], dtype = float)
with open(X_test_path) as f:
    next(f)
    X_test = np.array([line.strip(\'\\n\').split(\',\')[1:] for line in f], dtype = float)

#先把X_train 和 Y_train 组合成一个矩阵，构成一个数据集data，最后一列是真实标签，这一个操作用于方便后面的随机打乱操作
#(在Porbabilistic generative model中，随机打乱操作好像没什么用。。,后面也没有用随机打乱。)
Y_train=Y_train.reshape(-1,1)
data = np.concatenate((X_train,Y_train),axis=1)
print(data.shape)

# 标准化训练集, 这里没有划分验证集
train_set = data
train_set_mean = np.mean(data[:,0:-1] ,axis=0)
train_set_mean = train_set_mean.reshape(1,510)
print(train_set_mean.shape)
train_set_std  = np.std(data[:,0:-1], axis=0).reshape(1, -1)
print(train_set_std.shape)
train_set[:,0:-1]  = (data[:,0:-1]-train_set_mean)/(train_set_std+0.00000001) #避免方差为0。当方差为零时，说明这一列的数值是一样的，那就把这一列都变为0
print(train_set.shape)
print(train_set)

先看一下Porbabilistic generative model需要哪些参数

ppt中，告诉我们，计算出类别1的均值μ1和类别2的均值μ2，以及共享协方差Σ，就可以得到 w 和 b ，然后就可以计算 x 属于 类别1（在作业中，类别1对应标签0）的概率了，也就是P(C1|x)的概率。

看到这里，也清楚了Porbabilistic generative model是怎么一回事了。

Porbabilistic generative model 中的w和b是我们计算出来的，logistic regression 中的w和b是模型通过梯度下降自己学习出来的！ 这两个模型都是在找w和b！

接下来就开始计算μ1 和μ2 以及 Σ

先把类别1（标签为0）和类别2（标签为1）的数据分开

e0 = (train_set[:,-1]==0) #设置筛选条件，train_set最后一列真实标签为0的数据
e1 = (train_set[:,-1]==1) #设置筛选条件，train_set最后一列真实标签为1的数据
train_0 = train_set[e0]
train_0 = train_0[:,0:-1] #舍去标签值
train_1 = train_set[e1]
train_1 = train_1[:,0:-1] #舍去标签值

计算 μ1 和 μ2

注意：在ppt里，μ1 和 μ2 是列向量，但是这里的是行向量

mean_0 = np.mean(train_0, axis = 0)  #μ1
print(mean_0.shape)
mean_1 = np.mean(train_1, axis = 0)  #μ2
print(mean_1.shape)
print(mean_1)

接下来就是计算 协方差 Σ

首先要计算出类别1的协方差 Σ1 和类别2的协方差 Σ2

这里介绍一下协方差的计算方法。

一个就是直接用numpy的函数numpy.cov(x,rowvar=False))

另一个就是利用公式 cov =1/(N−1) * (X - μ).T· (X - μ)
(N表示样本数，T表示转置，后面的·表示矩阵相乘运算)

例如：

x = np.array([ [1, 2, 3,4],
          [3, 4, 4,5]])
c = np.zeros((4,4))
print(np.cov(x,rowvar=False))
me = np.mean(x,axis=0)
print(np.dot((x-me).T,x-me)/(2-1))）

[[2.  2.  1.  1. ]
 [2.  2.  1.  1. ]
 [1.  1.  0.5 0.5]
 [1.  1.  0.5 0.5]]
 
[[2.  2.  1.  1. ]
 [2.  2.  1.  1. ]
 [1.  1.  0.5 0.5]
 [1.  1.  0.5 0.5]]

计算协方差

# 分别计算类别0和类别1的协方差
cov_0 = np.zeros((len(mean_0), len(mean_0)))  #直接写成cov_0 = np.zeros((510, 510))
cov_1 = np.zeros((510,510))

cov_0 =  np.dot((train_0-mean_0).T,train_0-mean_0)/(len(train_0)-1)
cov_1 =  np.dot((train_1-mean_1).T,train_1-mean_1)/(len(train_1)-1)
print(cov_0.shape)

# 共享协方差 = 独立的协方差的加权求和
cov = (cov_0 * train_0.shape[0] + cov_1 * train_1.shape[0]) / (train_0.shape[0]+train_1.shape[0])
print(cov.shape)
print(cov)
print(np.cov(train_0,rowvar=False))

接下来就是计算 W 和 b 了

具体怎么计算，ppt里有公式，具体的推导过程可以自己参考ppt上的其他内容

在计算 w 和 b 时，需要计算 Σ 的逆矩阵

在numpy中可以直接用函数numpy.linalg.inv(a)求矩阵的逆，但这要求a是一个非奇异矩阵，我用这个函数计算Σ的逆矩阵，结果说Σ是一个奇异矩阵（秩不是满秩）。

下面的代码用numpy.linalg.svd()函数求解Σ的逆（Moore-Penrose 广义逆）。具体怎么求需要学习线性代数和矩阵论的一些知识，在这里直接用。

u, s, v = np.linalg.svd(cov, full_matrices=False)
inv_cov = np.matmul(v.T * 1 / s, u.T) #求Σ逆
print(inv_cov)

下面就可以计算 w 和 b 了

这里需要注意一下 ppt 中的μ是列向量，这里计算出来的是行向量，所以在计算时在转置的处理上有一点变化

# 计算w和b
w = np.dot(mean_0 - mean_1,inv_cov).reshape(-1,1)
b =  (-0.5) * np.dot(mean_0, np.dot(inv_cov, mean_0.T)) + 0.5 * np.dot(mean_1, np.dot(inv_cov, mean_1.T)) + np.log(float(train_0.shape[0]) / train_1.shape[0])
print(w.shape)

计算训练集上的准确率

Porbabilistic generative model 计算的是 x 属于 类别1（标签值为0）的概率，比如当p(c1|x) = 0.9时，此时x应该判别为类别1，打上标签0，

0.9四舍五入后是1，与标签0不对应，这时用 1- numpy.round(0) 就对应上了。

# 计算训练集上的准确率
x = data[:,0:-1]
x = x.reshape(len(data),-1)
y_label = data[:,-1] #获取标签值
y_label = y_label.reshape(-1,1)

f = np.matmul(x, w) + b
y_pred = 1- np.round(sigmoidFun(f)) #预测
acc = accuracy(y_pred,y_label)
print(\"训练准确率: \"+str(acc))

#标准化测试集
test_set  = (X_test-train_set_mean)/(train_set_std+0.00000001) #均使用训练集的均值和方差
print(test_set.shape)
print(test_set)

总结

在进行标准化操作时，这里是先划分训练集和验证集，然后在分别标准化，验证集标准化时，用训练集的均值和方差

Porbabilistic generative model 与 logistic regression 都需要找到w和b，Porbabilistic generative model中的w和b是 根据极大似然估计来确定 w 和 b ，就是模型已定，参数未知。 这里就是假设的分布是高斯分布，然后确定μ和Σ，最终确定w和b

logistic regression 模型通过梯度下降自己学习出来。