年收入判断

本项目是在飞桨深度学习学院提供的李宏毅-机器学习特训营课程。

李宏毅机器学习作业2-年收入判断

项目传送门

项目描述

二元分类是机器学习中最基础的问题之一，在这份教学中，你将学会如何实作一个线性二元分类器，来根据人们的个人资料，判断其年收入是否高于 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 这两个原始资料档则可以提供你一些额外的资讯。

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

特征格式

1. train.csv，test_no_label.csv。

• 基于文本的原始数据
• 去掉不必要的属性，平衡正负比例。

2. 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 " 。

项目要求

1. 请动手编写 gradient descent 实现 logistic regression
2. 请动手实现概率生成模型。
3. 单个代码块运行时长应低于五分钟。
4. 禁止使用任何开源的代码(例如，你在GitHub上找到的决策树的实现)。

数据准备

项目数据保存在：work/data/ 目录下。

!pip install --upgrade pandas

Looking in indexes: https://mirror.baidu.com/pypi/simple/
Collecting pandas
[?25l  Downloading https://mirror.baidu.com/pypi/packages/f3/d4/3fe3b5bf9886912b64ef040040aec356fa48825e5a829a84c2667afdf952/pandas-1.2.3-cp37-cp37m-manylinux1_x86_64.whl (9.9MB)
[K     |████████████████████████████████| 9.9MB 16.6MB/s eta 0:00:01
[?25hCollecting numpy>=1.16.5 (from pandas)
[?25l  Downloading https://mirror.baidu.com/pypi/packages/70/8a/064b4077e3d793f877e3b77aa64f56fa49a4d37236a53f78ee28be009a16/numpy-1.20.1-cp37-cp37m-manylinux2010_x86_64.whl (15.3MB)
[K     |████████████████████████████████| 15.3MB 9.1MB/s eta 0:00:011
[?25hRequirement already satisfied, skipping upgrade: pytz>=2017.3 in /opt/conda/envs/python35-paddle120-env/lib/python3.7/site-packages (from pandas) (2019.3)
Requirement already satisfied, skipping upgrade: python-dateutil>=2.7.3 in /opt/conda/envs/python35-paddle120-env/lib/python3.7/site-packages (from pandas) (2.8.0)
Requirement already satisfied, skipping upgrade: six>=1.5 in /opt/conda/envs/python35-paddle120-env/lib/python3.7/site-packages (from python-dateutil>=2.7.3->pandas) (1.15.0)
[31mERROR: blackhole 0.3.2 has requirement pandas<=1.1.5,>=0.24.0, but you'll have pandas 1.2.3 which is incompatible.[0m
[31mERROR: blackhole 0.3.2 has requirement xgboost==1.1.0, but you'll have xgboost 1.3.3 which is incompatible.[0m
Installing collected packages: numpy, pandas
  Found existing installation: numpy 1.16.4
    Uninstalling numpy-1.16.4:
      Successfully uninstalled numpy-1.16.4
  Found existing installation: pandas 1.1.5
    Uninstalling pandas-1.1.5:
      Successfully uninstalled pandas-1.1.5
Successfully installed numpy-1.20.1 pandas-1.2.3

# 导入相关库
import numpy as np
import matplotlib.pyplot as plt

1.Logistic regression

Logistic regression（逻辑斯蒂回归）是统计学中经典的分类方法。回归模型是根据输入特征可分为连续型和离散型。
当输入特征是连续型时，通过对训练数据的拟合，进而得到回归模型来预测房屋售价、产品销量等；而当输入特征是离散型时，用于解决分类任务。
例如判断一幅图片上的动物是一只猫还是一只狗例如，预测明天天气-阴，晴，雨。不仅适用二分类，还能解决多分类问题。Logistic regression
属于判别模型(Discriminative Model)通过学习得到条件概率分布P(y|x)，即在特征x出现的情况下标记y出现的概率。

Logistic regression的学习方法与算法:

Step1建立模型：

Step2采取的策略：

由此最大似然估计法等价与最小化交叉熵代价函数

Step3学习算法：通常采用的方法是梯度下降法和拟牛顿法，下面采用的是梯度下降法(Gradient Descent)

1.1 加载数据

直接导入已经处理好的数据X_train,Y_train,X_test

#添加文件路径
X_train_fpath = 'work/data/X_train'
Y_train_fpath = 'work/data/Y_train'
X_test_fpath = 'work/data/X_test'
predict_fpath = 'work/predict_{}.csv'	#用于测试集的预测输出

#加载数据
with open(X_train_fpath) as f:
    next(f)
    X_train = np.array([line.strip('\n').split(',')[1:] for line in f], dtype = float)
with open(Y_train_fpath) as f:
    next(f)
    Y_train = np.array([line.strip('\n').split(',')[1] for line in f], dtype = float)
with open(X_test_fpath) as f:
    next(f)
    X_test = np.array([line.strip('\n').split(',')[1:] for line in f], dtype = float)

# 查看数据
# X_train
X_train.shape,X_train.dtype,X_train

((54256, 510), dtype('float64'), array([[33.,  1.,  0., ..., 52.,  0.,  1.],
        [63.,  1.,  0., ..., 52.,  0.,  1.],
        [71.,  0.,  0., ...,  0.,  0.,  1.],
        ...,
        [16.,  0.,  0., ...,  8.,  1.,  0.],
        [48.,  1.,  0., ..., 52.,  0.,  1.],
        [48.,  0.,  0., ...,  0.,  0.,  1.]]))

# 查看数据
# Y_train
Y_train.shape,Y_train.dtype,Y_train

((54256,), dtype('float64'), array([1., 0., 0., ..., 0., 0., 0.]))

# 查看数据
# X_test
X_test.shape,X_test.dtype,X_test

((27622, 510), dtype('float64'), array([[37.,  1.,  0., ..., 52.,  0.,  1.],
        [48.,  1.,  0., ..., 52.,  0.,  1.],
        [68.,  0.,  0., ...,  0.,  1.,  0.],
        ...,
        [38.,  1.,  0., ..., 52.,  0.,  1.],
        [17.,  0.,  0., ..., 40.,  1.,  0.],
        [22.,  0.,  0., ..., 25.,  1.,  0.]]))

1.2 归一化处理

# 一般的方法，但是这种方法计算比较慢
# X_mean = np.mean(X_train,axis=0) # 每个特征的均值
# X_std = np.std(X_train,axis=0) # 每个特征的标准差
# n = X_train.shape[0]
# m = X_train.shape[1]
# X = np.full_like(X_train,fill_value=0)

# for i in range(n):
#     for j in range(m):
#         X[i,j] = (X_train[i,j] - X_mean[j]) / (X_std[j]+1e-6) # 归一化数据

# X

编写一个_normalize()函数对数据进行预处理：归一化，即每个数据特征的均值和标准差进行归一化

#归一化
def _normalize(X, train = True, specified_column = None, X_mean = None, X_std = None):
    # This function normalizes specific columns of X.
    # The mean and standard variance of training data will be reused when processing testing data.
    #
    # Arguments:
    #     X: data to be processed
    #     train: 'True' when processing training data, 'False' for testing data
    #     specific_column: indexes of the columns that will be normalized. If 'None', all columns
    #         will be normalized.
    #     X_mean: mean value of training data, used when train = 'False'
    #     X_std: standard deviation of training data, used when train = 'False'
    # Outputs:
    #     X: normalized data
    #     X_mean: computed mean value of training data
    #     X_std: computed standard deviation of training data

    if specified_column == None:
        #为每个数据添加索值
        specified_column = np.arange(X.shape[1])
    if train:
        #求取每个数据的平均值和标准差
        X_mean = np.mean(X[:, specified_column] ,0).reshape(1,-1)
        X_std  = np.std(X[:, specified_column], 0).reshape(1,-1)
    #归一化数据
    X[:,specified_column] = (X[:, specified_column] - X_mean) / (X_std + 1e-8)
     #返回归一化后的数据，均值，标准差
    return X, X_mean, X_std

1.3 划分数据集

利用_train_dev_split（）在train数据上分割出验证集，用来验证我们的模型

#分割训练集-验证集
def _train_dev_split(X, Y, dev_ratio = 0.25):
    # This function spilts data into training set and development set.
    train_size = int(len(X) * (1 - dev_ratio))
    return X[:train_size], Y[:train_size], X[train_size:], Y[train_size:]

现在我们来对X_train,X_test进行归一化处理
利用X_train归一化后的X_mean和X_std来处理X_test

X_train,X_mean,X_std = _normalize(X_train,train=True)
X_test,_,_ =  _normalize(X_test,train=False,specified_column=None,X_mean=X_mean,X_std=X_std)

采用dev_ratio = 0.1 的比例来设立训练-验证集

# 设置训练集-验证集
dev_ratio = 0.1
X_train, Y_train, X_dev, Y_dev = _train_dev_split(X_train, Y_train, dev_ratio = dev_ratio)
train_size = X_train.shape[0]
dev_size = X_dev.shape[0]
test_size = X_test.shape[0]
data_dim = X_train.shape[1]
print('Size of training set: {}'.format(train_size))
print('Size of development set: {}'.format(dev_size))
print('Size of testing set: {}'.format(test_size))
print('Dimension of data: {}'.format(data_dim))

Size of training set: 48830
Size of development set: 5426
Size of testing set: 27622
Dimension of data: 510

1.4 模型准备

训练集有48830个数据，验证集有5426个数据，测试集有27622个数据，数据维度为510

现在数据已经建立完毕，我们还需要定义一些函数：

• _shuffle：打乱数据顺序，类似于重新洗牌，进行分批次训练（即每次将一部分数据喂给模型进行训练，计算损失）
• _sigmoid：激活函数
• _f（X，w，b）：向前传播，计算激活值
• _predict(X, w, b)：预测
• _accuracy(Y_pred, Y_label)：计算准确度
• _cross_entropy_loss：交叉熵损失函数
• _gradient：计算梯度值，用于更新w和b
  代码实现如下：

#打乱数据顺序，重新为minibatch分配
def _shuffle(X, Y):
    # This function shuffles two equal-length list/array, X and Y, together.
    randomize = np.arange(len(X))
    np.random.shuffle(randomize)
    return (X[randomize], Y[randomize])

#sigmoid函数
def _sigmoid(z):
    # Sigmoid function can be used to calculate probability.
    # To avoid overflow, minimum/maximum output value is set.
    return np.clip(1 / (1.0 + np.exp(-z)), 1e-8, 1 - (1e-8))
#向前传播然后利用sigmoid激活函数计算激活值
def _f(X, w, b):
    # This is the logistic regression function, parameterized by w and b
    #
    # Arguements:
    #     X: input data, shape = [batch_size, data_dimension]
    #     w: weight vector, shape = [data_dimension, ]
    #     b: bias, scalar
    # Output:
    #     predicted probability of each row of X being positively labeled, shape = [batch_size, ]
    return _sigmoid(np.matmul(X, w) + b)

#预测
def _predict(X, w, b):
    # This function returns a truth value prediction for each row of X 
    # by rounding the result of logistic regression function.
    return np.round(_f(X, w, b)).astype(np.int)

#准确度
def _accuracy(Y_pred, Y_label):
    # This function calculates prediction accuracy
    acc = 1 - np.mean(np.abs(Y_pred - Y_label))
    return acc
#交叉熵损失函数
def _cross_entropy_loss(y_pred, Y_label):
    # This function computes the cross entropy.
    #
    # Arguements:
    #     y_pred: probabilistic predictions, float vector
    #     Y_label: ground truth labels, bool vector
    # Output:
    #     cross entropy, scalar
    cross_entropy = -np.dot(Y_label, np.log(y_pred)) - np.dot((1 - Y_label), np.log(1 - y_pred))
    return cross_entropy

#计算梯度值
def _gradient(X, Y_label, w, b):
    # This function computes the gradient of cross entropy loss with respect to weight w and bias b.
    y_pred = _f(X, w, b)
    pred_error = Y_label - y_pred
    w_grad = -np.sum(pred_error * X.T, 1)
    b_grad = -np.sum(pred_error)
    return w_grad, b_grad

1.5模型训练

# 至此，模型已经建立完成，我们开始训练
# 将w和b初始化为0
w = np.zeros((data_dim,)) 
b = np.zeros((1,))

# 设置其他超参数（迭代次数，分批次大小，学习率）  
max_iter = 100
batch_size = 128
learning_rate = 0.2

# 创建列表用来保存训练集和验证集的损失值和准确度
train_loss = []
dev_loss = []
train_acc = []
dev_acc = []

# 用来更新学习率
step = 1

# 训练
for epoch in range(max_iter):
    # 每个epoch都会重新洗牌
    X_train, Y_train = _shuffle(X_train, Y_train)
        
    # 分批次训练
    for idx in range(int(np.floor(train_size / batch_size))):
        X = X_train[idx*batch_size:(idx+1)*batch_size]
        Y = Y_train[idx*batch_size:(idx+1)*batch_size]

        # 计算梯度值
        w_grad, b_grad = _gradient(X, Y, w, b)
            
        # 更新参数w和b
        # 学习率随着迭代时间增加而减少
        w = w - learning_rate/np.sqrt(step) * w_grad
        b = b - learning_rate/np.sqrt(step) * b_grad

        step = step + 2
  	# 参数总共更新了max_iter × （train_size/batch_size）次     
    # 计算训练集的损失值和准确度
    y_train_pred = _f(X_train, w, b)
    Y_train_pred = np.round(y_train_pred)
    train_acc.append(_accuracy(Y_train_pred, Y_train))
    train_loss.append(_cross_entropy_loss(y_train_pred, Y_train) / train_size)
	# 计算验证集的损失值和准确度
    y_dev_pred = _f(X_dev, w, b)
    Y_dev_pred = np.round(y_dev_pred)
    dev_acc.append(_accuracy(Y_dev_pred, Y_dev))
    dev_loss.append(_cross_entropy_loss(y_dev_pred, Y_dev) / dev_size)

print('Training loss: {}'.format(train_loss[-1]))
print('Development loss: {}'.format(dev_loss[-1]))
print('Training accuracy: {}'.format(train_acc[-1]))
print('Development accuracy: {}'.format(dev_acc[-1]))

/opt/conda/envs/python35-paddle120-env/lib/python3.7/site-packages/ipykernel_launcher.py:12: RuntimeWarning: overflow encountered in exp
  if sys.path[0] == '':


Training loss: 0.2688210307575216
Development loss: 0.29058861468589503
Training accuracy: 0.8837599836166291
Development accuracy: 0.8730187983781792

1.6 训练过程可视化

绘制训练集和验证集的损失函数曲线，准确率曲线

# Loss Curve
plt.plot(train_loss)
plt.plot(dev_loss)
plt.title('Loss Curve1')
plt.ylabel('loss')
plt.xlabel('epoch')
plt.legend(['train', 'dev'])
plt.savefig('loss.png')
plt.show()

# Accuracy Curve
plt.plot(train_acc)
plt.plot(dev_acc)
plt.title('Accuracy Curve1')
plt.ylabel('acc')
plt.xlabel('epoch')
plt.legend(['train', 'dev'])
plt.savefig('acc.png')
plt.show()

1.6 输出预测结果

测试：将训练出来的参数应用到测试集上，得到预测结果,并保存在文件’predict_logistic.csv’中

# Predict testing labels
predictions = _predict(X_test, w, b)
with open(predict_fpath.format('logistic'), 'w') as f:
    f.write('id,label\n')
    for i, label in  enumerate(predictions):
        f.write('{},{}\n'.format(i, label))
 
# 打印一下数据前10项特征对应的权重
# ind = np.argsort(np.abs(w))[::-1]
# with open(X_test_fpath) as f:
#     content = f.readline().strip('\n').split(',')
# features = np.array(content)
# for i in ind[0:10]:
#     print(features[i], w[i])

/opt/conda/envs/python35-paddle120-env/lib/python3.7/site-packages/ipykernel_launcher.py:12: RuntimeWarning: overflow encountered in exp
  if sys.path[0] == '':
/opt/conda/envs/python35-paddle120-env/lib/python3.7/site-packages/ipykernel_launcher.py:29: DeprecationWarning: `np.int` is a deprecated alias for the builtin `int`. To silence this warning, use `int` by itself. Doing this will not modify any behavior and is safe. When replacing `np.int`, you may wish to use e.g. `np.int64` or `np.int32` to specify the precision. If you wish to review your current use, check the release note link for additional information.
Deprecated in NumPy 1.20; for more details and guidance: https://numpy.org/devdocs/release/1.20.0-notes.html#deprecations

2 Generative Model生成模型

Generative model方法跟Logistic regression方法类似，不同之处在于Generative model可以直接计算出w和b的最佳解，而Logistic regression是将w和b进行初始化，通过迭代训练来更新w和b

代码除了求解w和b地方不一样，其他地方类似。所谓生成模型（Generative Model）是模型表示了给定输入X产生输出Y的生成关系。由数据学习联合概率密度分布P(X,Y)，然后求出条件概率分布P(Y|X)作为预测的模型，即生成模型：P(Y|X)= P(X,Y)/ P(X)。

Generative Model的学习方法与算法:

Step1建立模型：Bayes理论

Step2采取策略：使用极大似然估计（使用样本中的数据分布来拟合数据的实际分布概率）得到先验概率和条件概率。

Step3学习算法：后验概率最大化

2.1 准备数据集

训练集与测试集的处理方法跟 logistic regression 一模一样，然而因为 generative model 有可解析的最佳解，因此不必使用到验证集（development set）。

# Parse csv files to numpy array
with open(X_train_fpath) as f:
    next(f)
    X_train = np.array([line.strip('\n').split(',')[1:] for line in f], dtype = float)
with open(Y_train_fpath) as f:
    next(f)
    Y_train = np.array([line.strip('\n').split(',')[1] for line in f], dtype = float)
with open(X_test_fpath) as f:
    next(f)
    X_test = np.array([line.strip('\n').split(',')[1:] for line in f], dtype = float)
 
# Normalize training and testing data
X_train, X_mean, X_std = _normalize(X_train, train = True)
X_test, _, _= _normalize(X_test, train = False, specified_column = None, X_mean = X_mean, X_std = X_std)

2.2 建立模型

具体的步骤如下：

• 在生成模型中，需要分别计算两个类别的每个特征的均值和标准差，
• 直接计算权重矩阵w和偏置向量b。

# 分别计算两个类别的每个特征的均值和标准差，
X_train_0 = np.array([x for x, y in zip(X_train, Y_train) if y == 0])
X_train_1 = np.array([x for x, y in zip(X_train, Y_train) if y == 1])
 
mean_0 = np.mean(X_train_0, axis = 0)
mean_1 = np.mean(X_train_1, axis = 0)  
 
# Compute in-class covariance
cov_0 = np.zeros((data_dim, data_dim))
cov_1 = np.zeros((data_dim, data_dim))
 
for x in X_train_0:
    cov_0 += np.dot(np.transpose([x - mean_0]), [x - mean_0]) / X_train_0.shape[0]
for x in X_train_1:
    cov_1 += np.dot(np.transpose([x - mean_1]), [x - mean_1]) / X_train_1.shape[0]
 
# Shared covariance is taken as a weighted average of individual in-class covariance.
cov = (cov_0 * X_train_0.shape[0] + cov_1 * X_train_1.shape[0]) / (X_train_0.shape[0] + X_train_1.shape[0])

权重矩阵w和偏置向量b可以直接计算出来：

# 计算权重矩阵w和偏置向量b
# Compute inverse of covariance matrix.
# Since covariance matrix may be nearly singular, np.linalg.inv() may give a large numerical error.
# Via SVD decomposition, one can get matrix inverse efficiently and accurately.
u, s, v = np.linalg.svd(cov, full_matrices=False)
inv = np.matmul(v.T * 1 / s, u.T)
 
# Directly compute weights and bias
w = np.dot(inv, mean_0 - mean_1)
b =  (-0.5) * np.dot(mean_0, np.dot(inv, mean_0)) + 0.5 * np.dot(mean_1, np.dot(inv, mean_1))\
    + np.log(float(X_train_0.shape[0]) / X_train_1.shape[0]) 
 
# Compute accuracy on training set
Y_train_pred = 1 - _predict(X_train, w, b)
print('Training accuracy: {}'.format(_accuracy(Y_train_pred, Y_train)))

Training accuracy: 0.8725486582129165


/opt/conda/envs/python35-paddle120-env/lib/python3.7/site-packages/ipykernel_launcher.py:29: DeprecationWarning: `np.int` is a deprecated alias for the builtin `int`. To silence this warning, use `int` by itself. Doing this will not modify any behavior and is safe. When replacing `np.int`, you may wish to use e.g. `np.int64` or `np.int32` to specify the precision. If you wish to review your current use, check the release note link for additional information.
Deprecated in NumPy 1.20; for more details and guidance: https://numpy.org/devdocs/release/1.20.0-notes.html#deprecations

2.3 输出预测结果

测试：将训练出来的参数应用到测试集上，得到预测结果,并保存在文件’predict_generative.csv’中

# Predict testing labels
predictions = 1 - _predict(X_test, w, b)
with open(predict_fpath.format('generative'), 'w') as f:
    f.write('id,label\n')
    for i, label in  enumerate(predictions):
        f.write('{},{}\n'.format(i, label))
 
# 打印一下数据前10项特征对应的权重
# ind = np.argsort(np.abs(w))[::-1]
# with open(X_test_fpath) as f:
#     content = f.readline().strip('\n').split(',')
# features = np.array(content)
# for i in ind[0:10]:
  enumerate(predictions):
        f.write('{},{}\n'.format(i, label))
 
# 打印一下数据前10项特征对应的权重
# ind = np.argsort(np.abs(w))[::-1]
# with open(X_test_fpath) as f:
#     content = f.readline().strip('\n').split(',')
# features = np.array(content)
# for i in ind[0:10]:
#     print(features[i], w[i])

/opt/conda/envs/python35-paddle120-env/lib/python3.7/site-packages/ipykernel_launcher.py:29: DeprecationWarning: `np.int` is a deprecated alias for the builtin `int`. To silence this warning, use `int` by itself. Doing this will not modify any behavior and is safe. When replacing `np.int`, you may wish to use e.g. `np.int64` or `np.int32` to specify the precision. If you wish to review your current use, check the release note link for additional information.
Deprecated in NumPy 1.20; for more details and guidance: https://numpy.org/devdocs/release/1.20.0-notes.html#deprecations

参考文献&文章&代码

(1) 李航.统计学习方法（第二版）
(2)https://blog.csdn.net/weixin_42447868/article/details/105294890
(3) https://qushuo.blog.csdn.net/article/details/106987801
(4)https://blog.csdn.net/weixin_39910711/article/details/89483662

作者介绍

大家好，我是黄波波。希望能和大家共进步，错误之处恳请指出！
百度AI Studio个人主页, 我在AI Studio上获得白银等级，点亮2个徽章，来互关呀~
交流qq:3207820044